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PROGRAMS FOR COMPUTING ABSCISSAS AND WEIGHTS FOR CLASSICAL 
AND NONCLASSIC AL GAUSSIAN QUADRATURE FORMULAS 

Robert N. Desmarais 
Langley Research Center 

SUMMARY 

Computer programs for computing Gaussian quadrature abscissas and weights are 
described. For the classical case the programs use Laguerre iteration to compute 
abscissas as zeros of orthogonal pol 3 raomials. The polynomials are evaluated from 
known recursion coefficients. The nonclassical case is handled similarly except that the 
recursion coefficients are computed by numerical integration. A sample problem, with 
input and output, is presented to illustrate the use of the programs. It computes the quad- 
rature abscissas and weights associated with the weight function (1 - x)^/^ In (1/x) over 
the interval (0,1) for quadrature orders from 16 to 96 in increments of 8. 

INTRODUCTION 

This paper describes a set of computer programs for computing the abscissas and 
weights of Gaussian quadrature formulas. The programs permit the calculation of both 
classical and nonclassical abscissas and weights. For the classical case the programs 
are complete in the sense that the user need only specify the order, interval, and constants 
appearing in the weight function, and the programs will do the rest. For the nonclassical 
case, the user has to set up the quadrature scheme to be used for computing the recursion 
coefficients of the orthogonal system associated with the Gaussian quadrature formula 
desired. 

The programs described herein were developed to generate quadrature formulas for 
use in computing unsteady aerodynamic forces. For example, the kernel of the integral 
equation relating lift to downwash in unsteady subsonic flow can be computed, in part, 
using a Laguerre -Gauss quadrature; integrals of the pressure distribution can be com- 
puted using Jacobi-Gauss quadrature; and integrals of pressure distributions over a 
control -surface hinge can be evaluated using a nonclassical Gaussian quadrature with the 
weight function In (1/x). 

For all of the quadrature formulas mentioned in the preceding paragraph except the 
one with the logarithmic weight the existing published tables of abscissas and weights are 



more than adequate. Even if tables are available, however, for classical weight functions 
it is more economical to generate the abscissas and weights from a computer program 
than it is to keypunch and verify a set of tabulated values. This is not true for nonclassi- 
cal weight functions. However, for these weight functions the existing tables are not 
adequate. 


SYMBOLS 


a,b limits of integration 

Cn>bn coefficients of recursion formula 

correction term for a singularity not appearing in the weight function 

f(x) function to be integrated 

fL(x) truncated Taylor’s expansion of f(x) about Xs 

g(x) either Pjj2(x) or xp^^2(x) 

hn integral defined by equation (14) 

hn integral defined by equation (15) 

I a definite integral to be approximated 

k index of summation in a quadrature formula 

L degree of truncated Taylor’s expansion 

L(x) Lagrange’s interpolation function 

I an index of summation 

7t(x) 

(X - Xk)7T’ 
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M 

m 

N 

n 

p(x) 

q(x) 

r(x) 

Sk 

Sg 

U(Q!,m) 

Wk 

X 

Xk 

Xs 

a,^ 

6nm 

y 

7t(x) 


number of terms in quadrature sum used to approximate hn or hn 

an index of summjition or degree of orthogonal pol 5 momial 
order of a Gaussian quadrature formula 

order of a Gaussian quadrature formula or degree of orthogonal polynomial 
orthogonal polynomials 

quotient when f(x) is divided by ir(x) 
remainder when f(x) is divided by ir(x) 

kth moment 

the fth shifted moment of p(x) 

= y + P{ot + m) 
a quadrature weight 

variable of integration 
a quadrature abscissa 

abscissa of a singularity of p(x) 

exponents appearing in p(x) or t(x) 

Kronecker delta (0 if n ^ m; 1 if n = m) 

Euler’s constant, 0,57721 . , . 

= (x - xi)(x - X 2 ) . . • (x - Xjj) where xjj for k = 1 to n are quadrature 
abscissas 
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p(x) weight function 

a(x) a factor of p(x) that is singular at x = Xg 

t(x) a classical or almost classical factor of p(x) 

i//(Z) digamma function, r’(Z)/r(Z) 

Prime denotes first derivative. 

Double prime denotes second derivative. 

Superscript within parentheses indicates a specific derivative. 

PROBLEM DESCRIPTION 

The most efficient Way to evaluate integrals with integrable singularities is to use 
a Gaussian quadrature formula. In this procedure the singular part of the integrand is 
factored out and processed analytically. Also the abscissas at which the remaining fac- 
tors of the integrand are to be evaluated are chosen so as to maximize the degree of the 
polynomial for which the procedure is exact. That is, a Gaussian quadrature formula is 
a numerical integration rule of the form 

b ” 

I = y p(x) f(x) dx « ^ Wk f(Xk) (1) 

^ k=l 


where p(x), Xk, and Wk are subject to the three following restrictions: 

1. The weight function p(x) does not change sign in (a,b). 

2. All integrals 

Sk = r p(x) x^ dx 


( 2 ) 


exist whenever k is a positive integer or zero. 

3. The quadrature abscissas xk and quadrature weights Wk are computed so that 
equation (1) is exact whenever f(x) is a polynomial of degree (2n - 1) or less. 

Given a weight function p(x) that satisfies conditions 1 and 2 above, the problem is 
to find the abscissas and weights that satisfy condition 3. In this paper a numerical pro- 
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cedure for computing these abscissas and weights is described along with a set of com- 
puter sul 5 )rograms written to facilitate implementation of the procedure. The procedure 
is based upon several well-known properties of orthogonal pol 3 momials. (See sec- 
tions 10.3 and 10.4 of ref. 1, for example.) These properties are: 

1. Gaussian quadrature abscissas xjj are the zeros of polynomials which are 
orthogonal with respect to the weight function p(x). 

2. The associated weights w^ can be computed from the abscissas and the orthog- 
onal polynomials. 

3. Consecutive orthogonal polynomials are connected by a three-term recursion 
relation. 

4. The recursion coefficients needed to evaluate a polynomial of any degree can be 
computed from integrals of products of lower degree orthogonal polynomials. 

A flow chart of the procedure is presented in figure 1. 


(start) 


COMPUTE "SIMPLER" 
ABSCISSAS ♦ WEIGHTS 
FROM KNOWN COEFF ICIENTS 


I SET LEADING RECURSION 
[COEFFICIENTS TOJ^; 

: 1 

EVALUATE INTEGRAND 
USING TABULATED RECURSION 
COEFFICIENTS AND PERFORM 
"SIMPLER” GAUSSIAN QUADRATURES 

I 

COMPUTE NEW RECURSION 
COEFFICIENTS FROM RATIOS OF 
INTEGRALS COMPUTED 
ABOVE AND STORE 
IN TABLE 



SOLVE POLYNOMIAL EQUATION 
FOR ABSCISSAS AND 
COMPUTE WEIGHTS THEREFROM 



Figure 1.- Flow chart of a method for computing Gaussian 
abscissas and weights using numerical integration to 
generate recursion coefficients of orthogonal polynomials. 
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Inspection of the flow chart shows that the recursion coefficients are computed by 
numerical integration and then a polynomial root finder is used to compute the abscissas 
and weights. 

The "simpler" Gaussian qxiadrature mentioned in the fourth box of the flow chart 
is a closely related quadrature (i.e., its weight function has some of the same singular- 
ities) for which the abscissas and weights are already known. Usually it is a classical 
Gaussian quadrature. It is used to remove as many as possible of the singularities of the 
weight function, and Taylor’s theorem is used to reduce the effect of those that remain. 
Because all singularities of an arbitrary weight function p(x) are either removed by 
incorporation into the "simpler" quadrature, or else have their effect reduced, this pro- 
cedure is very accurate. Because Taylor's theorem is used as a method of last resort 
to reduce the effect of singularities the procedure will be called the "Taylor's theorem 
method" in the rest of this paper. The mathematical details of Taylor's theorem method 
are described in the section entitled "Numerical Evaluation of hn and hn." 

Relation of Orthogonal Polynomials to Gaussian Quadrature 

Suppose that in equation (1) we have a set of abscissas for k = 1, 2, . . ., n 
chosen arbitrarily. Let 

7t(x) = (X - xi)(x - X2) . . . (x - Xn) (4) 


Then it can be seen that 


^k(x) = 


^(x) 

7T’(Xk)(X - Xk) 


(5) 


is equal to 1 if x = xjj and is equal to zero if x = x^, X2, . , or x^ (excluding x^). 
This means that the polsmomial 


L(x) = ^ f(xk) ^k(x) 
k=l 


( 6 ) 


coincides with f(x) at x = x^. The function L(x) is called the Lagrange interpolation 
polynomial. If L(x) is substituted for f(x) in equation (1) an expression for wjj is 
obtained 
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( 7 ) 



This expression is valid whether the xjj are Gaussian abscissas or are chosen 
arbitrarily. 

Now suppose that f(x) is a polynomial of degree (2n - 1). If f(x) is divided by 
7t(x), a quotient q(x) and a remainder r(x) both of degree (n - 1) are obtained. Thus 

f(x) = q(x) 77(x) + r(x) (8) 

If this is inserted into equation (1) the result is 


^b 

\ p(x) q(x) 7t(x) dx + \ p(x) r(x) dx 
<-'a »^a 



Wk q(xk) 7^(xk) + ^ wkr(xk) 

k=l 


(9) 


The first sum above is zero because 7r(Xk) =0 by equation (4). The approximation is 
exact for f(x), an arbitrary polynomial of degree (2n - 1), only if the first integral is also 
zero. This will occur if 7t(x) is orthogonal to all polynomials q(x) of degree (n - 1) or 
less. 

The abscissas and weights of a Gaussian quadrature formula can be obtained by 
constructing the sequence of polsmomials Pn(x) such that 

pb 

\ p(x) p^(x) Pj^(x) dx = 6nmhn (10) 


and then computing the zeros xjj of Pn(x) and computing 


pb p(x) p (x) 

Wu = \ dx 

-’a Pn(xk)(x - Xk) 


( 11 ) 


Recursion Formulas and Christoffel-Darboux Identity 

A set of polynomials Pn(x), n = 0, 1, 2, 3, . . ., for which equation (10) is true is 
called an orthogonal system with respect to the weight function p(x) over the interval 
(a,b). Equation (10) itself is not sufficient to define the polynomials Pn(x) uniquely. If 
hn is not specified p„ may have an arbitrary factor and if hn is specified the sign of 
Pjj is ambiguous. By specifying, in addition to equation (10), the coefficient kn of the 
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highest power term in Pn(x), p^ can be defined unambiguously. This is called the 
standardization of the system of polynomials. Except for the classical polynomials the 
standardization adopted herein is kn = 1. For this standardization the three -term 
recursion relation for a system of orthogonal polynomials is 


Pq(x) = 1 


Pl(x) = x - bi 




( 12 ) 


Pn(x) = (x - bn) Pn-lW - Cn Pn-2(^) 


where 


bn 





Cn ~ 


hn-l 

hn-2^ 


(13) 


and where 


hn = \ P(x) Pn^(x) dx (14) 

^a 

pb 

hn = \ P(x) X Pn^(x) dx (15) 

‘^a 

The expressions for bn and Cn are obtained by multiplying equations (12) through by 
P(x) Pn-l(x) and p(x) Pn_ 2 (x), respectively, and integrating from a to b. The same 
procedure using p(x) Pm(x), m = 0, 1, . . ., n - 3, furnishes an inductive proof that the 
sequence of polynomials generated by equations (12) is orthogonal. 

The integrals hn and hn are computed numerically and used to generate a table 
of recursion coefficients bjj,cij. If equations (12) are differentiated, recursion formulas 
for p^ and are obtained as follows: 


Pq(x) = 0 


(16a) 
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Pi(x) = 1 


(16b) 


p|j(x) = (x - bn) Pn_i(x) - Cn Pn_2(^) + 


ip^”>W = o 
ip<i”'w = o 






— pVj /(X) - (X - bn) - - Cn — 


m: 


(m - 1): 


(16c) 


(17) 


The quadrature abscissas x^ are obtained by solving the polynomial equation 
Pjj(x) = 0 using Laguerre's iteration formula as described in appendix A, Laguerre's 
iteration formula requires values of Pn(x), Pn(x), and p'n(x) and these are furnished 
by equations (12), (16), and (17), 

After the abscissas have been computed, the weights are computed from 


Wk 


hn-1 

PU^'k) Pn-l(^k) 


This is obtained from equation (11) as follows. Let 


Fn(x,y) 


Pn(x) Pn-l(y) - Pn-l(^) Pn(y) 
(x - y)hn-l 


(18) 


(19) 


If equations (12) are used to eliminate Pj^(x) and Pj^(y), one obtains 


Fn(x,y) = Pn-l(x) Pn-l^y) + Fn_i(x,y) 
%-l 


( 20 ) 


Repeating the process (n - 1) times leads to the sum 


n-1 


Fn(x,y)= ^ iPm(x)Pm(y) 


( 21 ) 


m=0 
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Equations (19) and (21) are a form of the Christoff el -Darboux identity for orthogonal 
polynomials (see eq. 10.3 (10) of ref. 1 for the more usual form). If y is replaced by 
Xk, a zero of Pn(x), the result is 


Pn(x) 
X -Xk 


^n-1 

Pn-l(^k) 


m=0 


( 22 ) 


Substituting this into equation (11) gives equation (18). Observe that all terms of the sum 
except the m = 0 term integrate to zero because of the orthogonality of the polynomials 

PmW- 


Numerical Evaluation of hn and hn 

The quantities hn and h^ used to compute the recursion coefficients are 
obtained by a combination of numerical and closed-form integration. The first step is to 
evaluate the first few moments of the weight function anal 3 dically. That is, 



(23) 


is computed for jf = 0, 1, 2, . . ., L. The integer L is the order of the singularity 
extraction (to be described later) and is usually 2 or 3. Then the first few pairs of 
recursion coefficients bk and ck are computed from these moments using equa- 
tions (13), (14), and (15) 


ho = Sq h{j = 


bl 



hi — S2 - 2bjSj + bj^^So 


‘^2 


= !ll 
ho 


> 

hj “ Sg - 2bjS2 “t" 


(24) 


h2 



J 
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This furnishes the entries up to k = L - 1 in a table of recursion coefficients 
bjjjCjj. Each time a pair of integrals hij,hJj is computed a new pair of entries 
b^+i = h{j/hjj and = hjjy/hjj,^ is added to this table. As a consequence, each time 

an hjj or hj^ integral is to be computed, all the recursion coefficients needed to eval- 
uate the integrand are available. 

Either h^ or h^ is obtained by evaluating the integral 


pb 

I = \ p(x) g(x) dx . (25) 

^a 

numerically where g(x) = p^j^(x) if h^ is being evaluated and g(x) = x Pj^^(x) if h^ 
is being evaluated. Since the objective of this report is to describe a procedure for com- 
puting quadrature abscissas and weights for a weight function p(x) that has a certain set 
of singular points, particular attention is paid to these singular points when equation (25) 
is integrated. Let p(x) be factored into two functions 

p(x) = ct(x) t(x) (26) 

where t is chosen so that (1) t(x) contains as many singularities of p(x) as possible, 
and (2) t(x) is a weight function over (a,b) whose Gaussian quadrature abscissas and 
weights are relatively easy to compute. Occasionally it will be necessary to let t(x) = 1. 

After p(x) has been factored the Gaussian quadrature abscissas x^ and weights 
W{j associated with the numerical integration 



k=l 


(27) 


are computed. The quadrature order N is usually taken to be about five times the max- 
imum order desired for the weight function p(x). 

Reduction of effect of singularities in o(x) . - If a(x) has no singularities then hn 
and hn are evaluated using equation (27) with f(x) = a(x) g(x). If a(x) has one or 
more singularities (the usual case) the strength of these singularities can be reduced by 
integrating the leading terms of Taylor's expansion about these singularities in closed 
form. Before describing the method in detail an example will be discussed showing the 
motivation for the method. 
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Consider the four integrals 


and 


4 = 



dx = 7T = 3.141 592 65 


- dx = f = 1.570 796 33 

Ck 

a 

- x2)^/^ dx = ^ = 1.178 097 25 


- x2)^/^ = If = 


If these integrals are approximated by eight -point Legendre -Gauss quadratures (i.e., 
eq. (1) with n = 8 and p(x) = 1) then 


G8(Ii) = 2.936 842 06 
G8(l2) = 1.572 152 22 
G8 (Is) = 1.178 033 49 


for a -6. 5 -percent error 
for a 0,087-percent error 
for a -0,005-percent error 


and 


^8(14) ~ 0.981 755 58 for a 0.0008 -percent error 


The reason for the difference in accuracy is the order of the singularities at ±1. The 
integrand of Ij is infinite at ±1 while I4 merely has an infinite third derivative. This 
is so because for each singularity Xs in the integrand of Ij the integrand of I4 has 
(x - Xg)^ as a factor. 

Now consider the integral I of equation (25), It can be written 


12 


pb 

I = \ t(x) a(x) g(x) dx 


( 28 ) 



It is assumed for now that a(x) has only one singularity at x = Xs, The detrimental 
effect of this singularity on the numerical quadrature of equation (27) is greatly reduced 
if g(x) is replaced by something that has a power of (x - Xs) as a factor. By Taylor's 
theorem 


g(x) = gL(x) + Rl(x) 


(29) 


where 


±i 

gL(x) = 


(X - Xg)^ 


g 


«){xs) 


jf=0 


(30) 


is the leading part of the Taylor's expansion of g(x) about Xg, the singularity of ct(x) 
and 


Rl(x) = g(x) - gL(x) 


(31) 


is the remainder. Note that Rl(x) has (x - Xg)^"*"^ as a factor. Equation (28) can be 
expressed 


pb pb 

I = J p(x) gL(x) dx + J r(x) a(x) Rl(x) dx 

Si Si 


(32) 


The first integral can be evaluated analytically 


pb 

II = j p(x) g^(x) 


dx 


i-i 

-I 




£=0 


i\ 




(33) 


where the shifted moments Sj^ are obtained from the moments about the origin that 
were previously computed (see eq. (23)) 


2 (-l)“(A)S£-mXs“ 

m=0 


(34) 
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The second integral in equation (32) is evaluated numerically using equation (27) 





( 35 ) 


Equations (33) and (35) can be combined to give 

N L 

I * ^ Wk a(xk) g(Xk) + ^ g(^)(xs) 

k=l H=0 


(36) 


where 


N 

Cl = S£ - ^ Wk a(Xk) (Xk - xs)- 


k=l 


(37) 


is independent of g(x). This is potentially a very accurate way to approximate an inte- 
gral, but care has to be exercised in choosing L. If L is larger than the degree of 
g(x) the quadrature is exact; that is, there is no truncation error. 1 In this case it is 
merely a scheme for expressing I as a linear combination of the moments and this 
method is known to be ill-conditioned (i.e., uncontrolled growth of round-off error^). 

Each correction coefficient Cj^ is the difference between the shifted moment 

pb 

Sjg = \ p(x) (x - Xg) ^ dx 
•^a 

and a Gaussian approximation to the shifted moment. As either 1 or N becomes 
large, will become very small. Since Ci is the difference of two quantities that 
are not approaching zero, eventually for some 1, Ci will have no significant figures 
at all. Thus, L should be chosen so that this does not happen. When this method is 
used, a plot of log jCil against 1 gives points that lie approximately on a straight line 
when Cji has significant figures and lie above the line when Cji has no significant fig- 
ures. Typical values for L are 2 or 3. 


^Truncation error is the error resulting from using insufficient terms in a limiting 
process, such as a series summation or a quadrature, while round-off error is the error 
resulting from performing mathematical operations with finite -length computer words. 
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If Pn^(x) and x Pn^(x) are substituted for g(x) in equation (36), the following 
numerical integration formulas for hn and are obtained: 


N L 

hn “ ^ Wk cr(xk) Pn2(xk) + ^ Cj qj^(xg) 
k=l £=0 




N L L 

hh “ 2 Wk cr(xk) Xk Pn^(xk) + xs ^ qje(xs) + ^ q£_i(xs) 
k=l H=0 i=l 


(38) 


where 


q^(x) = 


1 d^ 
dx^ 



pj”^^(x) p(£-m)(x) 
ml (fi - m)I 


(39) 


In equations (28) to (35) it was assumed that cr(x) had only a single singular point xg. 

If a(s) has a second singularity then Rl(x) = g(x) - g£j(x) must be expressed as a 
truncated Taylor’s series plus remainder expanded about the second singularity. 

Simpler quadratures .- An essential step in the evaluation of hn and hh for a 
particular weight function p(x) is the computation of the abscissas and weights associ- 
ated with a simpler weight function t(x) where t(x) is a factor of p(x). This so- 
called "simpler quadrature" that is indicated in equation (27) must have abscissas and 
weights that are computable without the necessity of evaluating their orthogonal poly- 
nomial recursion coefficients by nonexact numerical quadratures. Simpler quadratures 
can be classified as either classical or almost classical Gaussian quadratures depending 
upon the nature of their associated orthogonal polynomials. 

The classical orthogonal polynomials are the ssrtems of orthogonal polynomials that 
can be generated from a generalized Rodrigues' formula 

where is a constant and Q(x) is a pol 3 momial in x that is independent of n. It 
can be shown (see section 10.7 of ref. 1, for example) that the only zeros of Q(x) are the 
limits of integration in equation (1) and hence Q(x) must be of degree 2, 1, or 0. The 
associated orthogonal polynomials are, except possibly for a linear change of scale, the 
classical Jacobi, generalized Laguerre, and Hermite polynomials, respectively. If equa- 
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tion (40) is substituted for one of the Pn(x) factors in equations (14) or (15) the integral 
can be evalxiated in closed form by integration by parts. This furnishes closed form 
expressions for the recursion coefficients (see section 22.7 of ref. 2, for example). The 
classical orthogonal polsmomials also each satisfy a second-order differential equation 
and a first-order differential relation with respect to degree (see sections 22.6 and 22.8 
of ref. 2). This makes it possible to compute Pn(x) and pj{(x) directly from the 
recursion formula for Pn(x). Separate recursion relations such as equations (16) 
and (17) are not needed for classical polynomials. The fact that p^ and Pn are 
obtained free for the classical polynomials is the motivation for using Laguerre iteration 
for computing quadrature abscissas instead of the QR algorithm as in references 3, 4, 
and 5. 

The almost classical orthogonal polynomials are those for which equations (14) 
and (15) can be evaluated as exact classical Gaxissian quadratures. A Gaussian 
quadrature 



(41) 


is exact (i.e., no truncation error) if f(x) is a polynomial of degree (2n - 1) or less. 
Inspection of equations (14) and (15) shows that h^ and h^ can be evaluated by exact 
quadratures if p(x) can be expressed as a classical weight function times a polynomial 
or if (a,b) can be partitioned into a set of abuttii^ intervals such that within each interval 
p(x) can be expressed as the product of the classical weight function for that interval 
times a polynomial. Sometimes it is possible to replace p(x) by its integral definition 
and perform a multidimensional exact quadrature. For example, consider the weight 
function ln(l/x) over interval (0,1). Then 



In this case p(x) has the integral definition 



(42) 


(43) 
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so 




The substitution x = uv gives 


pi pi 

hh = \ \ uv p„2(uv) du dv 

■Jq ^0 


( 44 ) 


(45) 


This can be evaluated as an exact Gaussian quadrature 
n+1 n+1 

hn= 2 WkXk ^ WmXm Pn^(xkXm) (46) 

k=l m=l 

where xk,Xm,wk,Wni are the abscissas and weights associated with the (n + 1) -point 
classical Gaussian quadrature with weight function p(x) = 1 over (0,1). Several of 
these multidimensional exact quadratures are considered in reference 6. 

The purpose of the simpler quadrature is to minimize the number of singularities 
of p(x) = a(x)r(x) that have to be treated by Taylor's theorem. The simpler quadrature 
is the quadrature associated with the factor t(x) as a weight function and should be 
either a classical or almost classical Gaussian quadrature. 

PROGRAM ORGANIZATION 

The FORTRAN computer program to implement the Taylor's theorem method of 
calculating Gaussian quadrature abscissas and weights consists of a user-written calling 
prc^ram and a furnished subprogram package. The user-written calling program handles 
input/output (I/O) and that portion of the computing task that is peculiar to the weight 
function being processed. The subprogram package handles the portion of the program- 
ing task that is common to all weight functions. The complexity of the user -written 
calling program depends upon the nature of the weight function and the I/O services 
desired. For a classical weight function it could be as simple as a call to subroutine 
CGAUSS followed by a print statement whereas for a nonclassical weight function with 
several singularities in p(x) it could be very complicated. 

The subprogram package consists of four subprograms that are called by the user 
plus eight other subprograms. The former are: 
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CGAUSS a subroutine to compute classical Gaussian abscissas and weights for an 
arbitrary interval and arbitrary weight function exponents 

PNDER a subroutine to compute Pn”^\x) for m = 0,M from the recursion 

m. 

formulas (12), (16), and (17) 

PNFUN a function to compute Pn(x). This is an abridged version of PNDER that 
executes faster when only p^^^(x) is needed 

NGAUSS a subroutine to compute abscissas and weights for a nonclassical Gaussian 
quadrature 

For computing nonclassical Gaussian quadratures, the user-written calling program 
has two tasks to perform. The first and more difficult task is to compute the recursion 
coefficients b^ and Cn required by NGAUSS. ^ecifically this requires FORTRAN 
instructions to: 

1. Implement equations (23), (24), and (34). 

2. Compute the abscissas and weights required by equation (27) (a call to CGAUSS). 

3. Implement equation (37). . 

4. Compute as many bn,Cn as needed from equations (38), (39), and (13). Note 
that PNFUN and PNDER are used when implementing equations (38) and (39). 

The second and simpler task, to compute the desired nonclassical Gaussian abscissas 
and weights, is accomplished by simply a call to NGAUSS. The user-written calling pro- 
gram also manages two labeled COMMON blocks BOFN and COFN that contain the recur- 
sion coefficients bn and Cn that are given by equations (13). These coefficients are 
used by PNDER, PNFUN, and NGAUSS. NGAUSS also requires hg (see eq. (14)) and 
since cj is not used, hg is passed as the first word of COFN. The easiest way to 
describe in detail how a calling program is written is to explain on a step-by-step basis 
how the program was written for the sample problem. 

A Sample Problem Illustrating How Nonclassical Abscissas 
and Weights are Computed 

The sample problem is a FORTRAN program, calling CGAUSS, PNFUN, PNDER, 
and NGAUSS, to compute nonclassical Gaussian quadrature abscissas and weights to 
evaluate 
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( 47 ) 




n 

(1 - x)^/2 In 1 f(x) dx » ^ Wk f(xk) 

k=l 

for values of n from 8 to 96 in increments of 8. The weight function (1 - x)^/^ In (1/x) 
occurs when integrating the chordwise pressure distribution over an aircraft control sur- 
face in subsonic flow. The coordinates have been rescaled so that x = 0 is the location 
of the hinge line and x = 1 is the control-surface trailing edge. 

One of the purposes of the user -written calling program (henceforth referred to as 
prc^ram SAMPLE) is to handle input/output. It is possible for SAMPLE to have no input 
(i.e., all parameters are built in). However, a program of this sort should be fairly gen- 
eral, hence, as many parameters as possible are read in as input. 

In the problem statement min(n), max(n), and An were specified as 8, 96, and 8, 
respectively. These are read in as Nl, N2, and N3. Similarly, the L appearing in 
equation (36) is read in as LMAX, the number of correction coefficients computed 

and printed, and also as LX, the number of correction terms used in equation (39). The 
value of N appearing in equation (27) is read in as NC. The total number of pairs of 
recursion coefficients to be computed is read in as NMAX. 

Instead of writing the program to process the actual weight function specified in 
equation (47), it is written to process the more general function 

p(x) = (1 - x)"-l Inl (48) 

and a is read in using the program variable name ALPHA. A summary of the program 
input is as follows: 

IFLAG if 0, compute bn and Cn; if 1, read bn and Cn 
ALPHA a 

UA y + xl/{a) 

LMAX number of correction coefficients computed 

LX number of correction coefficients used 

NC order of classical Gaussian quadrature 
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NMAX 


number of bn,Cn pairs to be computed or read 


N1,N2,N3 delimiters of loop that calls NGAUSS 

Before program SAMPLE can be written, e?q)ressions for S£, the moments of p(x), 
must be derived and p(x) must be factored into a(x) and t(x) 

Sjf = \ (1 - x)“"^ In i x^ dx (49) 

Jq X 

To evaluate this integral, let x = 1 - u and note that 


1 ^ r" dv 

1 - u Jq 1 " V 


Then 






m=0 


m=0 


where 


(50) 


(51) 


pi 1 ,,Q!+m 

= i i -v ~ 


(52) 


This can be expressed as a digamma function (see eq. 6.3.22 of ref. 2) 
U(o!,m) = y + i//(q! + 1 + m) 


(53) 


which can be written (eq. 6.3.6 of ref. 2) 


m 

U(o!,m) = y + xl/{a) + ^ 


n=0 


(54) 
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Equation (54) suggests that U(Q!,m) should be computed recursively 


U(q!,-1) = y + if/(a) 


U(a,m) = U(a,m-1) + — i — 
\ > / \ > / 0! + m 


( 55 ) 


Instead of having the program compute the quantity U(ar,-1) = y + i//(o!) is read 

in with the variable name of UA. This is done because the digamma function takes on a 
simpler form when its argument is an integer or half integer. Some values of U(o!,-l) 
are 

U(l/2,-l) = -2 In 2 = -1.386 294 361 119 890 618 834 464 243 
U(l,-1) = 0 


U(3/2,-l) = 2 - 2 In 2 = 0.613 705 638 880 109 381 165 535 757 
U(2,-l) = 1 

U(5/2,-l) = I - 2 In 2 = 1.208 372 305 546 776 047 832 202 424 
o 

U(3,-l) =1 




(56) 


The weight function p(x) factors into (1 - x)““^ and In (1/x). Either could be used 
for the simpler quadrature of equation (27) because (1 - x)®"^ is classical and In (1/x) 
is almost classical. The choice 


ct(x) = In i 


r(x) = (1 -x)“-l 


(57) 

(58) 


leads to a simpler program. 

All the information needed to write the program is now available. Figure 2 is a 
flow chart of the program. Figure 3 is an expanded version of the box of the figure 2 
flow chart that contains the bn and cn calculations for n > 3. The prc^ram was 
written from the flow chart for the CDC .6000 series FORTRAN "RUN” compiler. The 
prc^ram listing follows along with a listing of subroutine COMPARE that was used for 
accuracy estimation when ALPHA = 1. Usage description of the subroutine package is 
given in appendix B. 
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Figure 2.- Program flow chart. 
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Sample -Problem Program Listing 


PROGRAM SAMPLE(INPUT=l ,OUTPUT=l ,PUNCH=l 
+ ,TAPE5= INPUT, T A PE6=OUT PUT, TAPP 7= PUNCH) 

CnMMON//X( 2000) ,W( 2000) 

00U8LF X,W 

common/ B 0FN/B( 100) /COFN/CdOO) 

DOUBLE B,C 

DOUBLE S( 10) ,Pb( 10) ,PS( 101 tOUO) 

DOUBLE ALPHA,UA,M,HP,HSf ADDtFAC 

DOUBLE DLOGfPNFUN 

READ 101, IFLAG 

IF( IFLAG. NE.O) GO TO 12 

* IF IFLAG. NE.O THE RECURSION COEFFICIENTS ARE RFAO IN. 


♦ PART 1 - RECURSION COEFFICIENT CALCULATION. 

RPAD 102, ALPHA, UA 
READ 101, LMAX 

COMPUTE U(ALPHA,M-1) AND STORE IN SIM) USING FQN (55) 

S( 1 ) =UA+1. /ALPHA 

no 1 m= 2 ,lmax 

1 S( M) = S( M-1M*1./(M-1+AL PH A) 

no 2 m=i,lmax 

2 S(M) = S(M)/(M-H'ALPHA) 

SGM=(-I.)««(LMAX-1) 

DO 4 LL=2,LMAX 

COMPUTE S(U FROM EON (51). 

L=2+LMAX-LL 
FAC = U 
DO 3 M=2,L 

FAC=-FAC^( I .+L-M)/ (M-1 . ) 

3 S(L I^SID+FAC^^'SI dL-M) 

S(L)=SGN*S(L) 

4 SGN=-SGN 

COMPUTF H0,B1,C2,AND B2 FROM FQN (24). NOTE THAT HO IS STORED IN C(l) 
C( 1) = S( 1) 

B(1 ) = S(2)/C( 1) 

H = S(3)-B(l)«(2.<^S(2)“B(l)*S(l)) 

C (2)=M/C( 1 ) 

B(2)=(S(4)-B(1)*(2.*S(3)-B(1)*S(2)))/H 
CALL CGAUSS TO SET UP SIMPLER GAUSSIAN QUADRATURE 
READ 101, NC 

CALL CGAUSS(NC,X,W,0.0, 1 . D , A L PH A« 1 . D , 0. D ) 

COMPUTF C(L-l) FROM FQN (37) AND STORE IN S(L). 

DO 5 K=1,NC 

ADO=W (K ) *OLOG(X (< ) ) 

S(1)=S( D+ADO 
DO 5 L=2,LMAX 
ADD = A00’4^X( K) 

5 S (L )=S( D+ADD 
PRINT 201, ALPHA, NC 
M = 0 

no 6 L=1,LMAX 

CLA = ALOG10( ABS( SMGLISIL) ) ) ) 

PRINT 202, M,S(L),M,CLA 

6 M=M+1 

READ 101, LX,NMAX 

COMPUTE RECURSION COF FF IC lENTS- FOR N=3 THRU NMAX . 

DO 10 N’ = 3,NMAX 

CALL PNDER{N-1,LX,0.D,PD,PS ) 
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COMPUTE 0 FROM EQN (39). 

DO 7 L=l.LX 

0( L)=0. 

DO 7 M=1,L 

7 0(L)=0(L)+PO(M)«Pn(L+l-M) 

HS = H 

H = SU)*0(1) 

HP=0. 

COMPUTE TAYLOR’S SERIES CORRECTION PART OF H,HP INTEGRALS IN EON (38). 
DO 8 L=2.LX 
H=H+S(L)*0(L) 

8 HP=HP+S(L)*0(L-1) 

DO 9 K=1,NC 

COMPUTE GAUSSIAN PART OF H,HP INTEGRALS IN ='QN (38). 

ADD=-DLOG(X(K ) ) *W ( K ) *PNFUN ( X(K),N-1)**2 
H=H+ADO 

9 HP=HP+ADD*X(K) 

B(N)=HP/H 

10 C (M)=H/HS 

PUNCH 301. ALPHA, LX, NC.NMAX 
DO 11 N=1,NMAX 

11 PUNCH 302, N,B(N),C(N) 

GO TO 14 


* PART 2 - QUADRATURE ABSCISSA AND WEIGHT CALCULATION. 

12 READ 301, ALPHA, LX, NC.NMAX 
no 13 N=1,NMAX 

13 READ 302, N,6(N) ,C(N) 

* THIS PART OF THE PROGRAM COMPUTFS THE ABSCISSAS AND WEIGHTS. 

* IT IS WEIGHT FUNCTION INDEPENDENT. 

14 READ 101, N1,N2,N3 
IF(N2.LE.0) STOP 2 
IF(N3.LF.0I L3=l 

PUNCH 303, ALPHA, N1,N2,N3 
DC 15 N=N1,N2,N3 
PUNCH 304, N 
CALL NGAUSS (N,X,W) 

DO* 15 M=1,N 

15 PUNCH 302, M,X(M),W(M) 

I F( ALPHA.EQ. 1.. AND.NMAX.GE.16 ) CALL COMPARE 
STOP 1 

101 FORMAT! 314) 

102 FORMAT (2040.30) 

201 F0RMAT(//36H CORRECTION COEFFICIENTS FOR ALPHA =,038.28/ 

+ 9H AND NC =, 15// ) 

202 F0RMATI4H C(,I1,3H) = , D38. 28 , 8 XI 2HLPG1 0 ( ASS ( C ( , 1 1 , 5H ) ) ) =,F8.2) 

301 FORMAT(//80H RFCURSION COEFFICIENTS FOR THE WEIGHT FUNCTION LN(1./ 
+X)*( 1 .-X)** ( ALPHA-1. ) 0VER/22H THE INTERVAL ( 0 . , 1 . ) ./ / 6X7HAL PHA =, 
+P36.28/9X4HLX =,I4/9X4HNC = , I 5/ 7X6HLMA X =,I4//3X1HN, 

+22X4HBIN) ,34X4HC(N) /) 

302 FORMAT! 14,2038.28) 

303 F0PMAT(//80H ABSCISSAS AND WFIGHTS FOR THE. WEIGHT FUNCTION LN(1./ 
+X)*( l.-X)**( ALPHA-1.) 0VER/34H THF INTERVAL (0.,1.) FOR ALPHA =, 
•FD3&.28,4H AN0/9H FOR N =13, 9H THRU N =I3,17H IN INCREMENTS OF, 13) 

■ 304 FORMAT! ///38X3HN = , I 3 //3X IHM, 22 X4HX ( M ) , 34X4HW ( M ) / ) 

END PROGRAM SAMPLE 
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SJB&njT IN*^ CGAUS'^C^ »X,W# A, Bt ALf’H^r3ETA I 
C'lMPUTHS ARS:rSSAS ANO a'B IGHT5 f=rjR CLASSICAL GUASSIAN .^JUADPATURE 
nnjBL R X ( I ) , W ( 1 ) , a, A, ALf’ HA , BETA 

:»I3JBL= co,ci rr^,os'^KT 

CASE I . ( JAC3BI ) 

1 IR ( BETA . L " 1. ) GO TU 3 

I r ( ALPH A .LF. -1 . ) GO r 3 5 
CALI J3AU5S(N,X,W, ALPHA, ItTA) 

CO=.5*( B+A) 
r 1= .s*( p-A ) 

C2 = r ( ALPHA-I-BETA+ l . ) 

D3 2 1= I tN 
W(I l=C2>^W( I) 

2 X( I ) = CO+CL^X( T ) 

5 c T UP N 

CASE 2. (LAGJFRRE) 

3 IF ( ALPHA. LE.-l . ) GJ TJ B 
CALL LGAUSS( \*f X,W, alpha ) 

: I =1. /3 

C? = Cl*<' ( AL PHA+l . ) «r)FX-> ( - 3>!^A ) 

30 4 1= l,N 

V n )=C2^w( n 

4 X( I)=A + Cl*X( I ) 

PFTUPN 

CASE 3 . ( HFPH I TC ) 

5 *MU=N/2 
N‘^ = NU+'4U 
Nl = H4»l 

IF(N.Gr.NNJ) Gn n 8 
"'ALL IGAUSS(NUfX,W,-,3 0) 

00 6 HJ=l,NJ 
X ( NIHHJ ) =os.;rt ( X ( ‘HU ) ) 

5 W ( Nl.i + MU ) = . 5*W( MU ) 

30 7 MU=ltMO 
X { MU) ^-X ( Nl-MJ ) 

7 MU) =W( Nl-MIJ) 

GO TO 13 

8 IF(NU.EQ.O) GO TO 12 
CALL LGAUSS( MU,X,W, .50) 

*viUl=NJ+ I 

3C 9 MU=l ,N'i 
X(NUH-MU) = DSQPT(X{MU) ) 

9 w ( NJU14-MIJ) =.54^W( MU) /X( MU) 
on 10 MJ=1,MU 

X( MU) =-X(Nl-MU) 

10 W(MU) = W(NI-M1J) 

X( MUl ):=0. 

W( NUl) = I . 7 72 453 8509 05 5 1502 7 29816 748333 
Cl = 2, 

C2 = 3. 

DO 11 MJ=l,MU 
W(NIU1)=W(NU1)«CI/C2 
Cl = CH-2 . 

U C2=:C2<'2. 

GO TO 13 

12 X=0. 

W = l.772 453 85 09055 15 0272981 6 74833D 

13 C1=1./3S0RT{ A) 

C0=.5*B/A 

C2=DEXP(A=«'C0*C0) 

00 14 M=l,N 

X(M ) = C1*X IM)-C0 

14 w(M)=C2^W(M) 

RETURN 
END 


SUBROUTINE CGAUSS 



SU=i^3UTINF JGAUSS(M ,X,W,<^L PH A, BETA) 

C")MPUTFS JACOB l-GAIJSS ABSCISSAS WEIGHTS F3B AL PHA , BE T A .GT 1. 
n:iJ3tE X( l ),W( DtALPHAfBETA 

OOUBLF XC ,P, PI f CON, A,H rC,0, Al,8l , A2, 82f AB2, StDSQRT 
CALL JWTCON(CnN,N» ALP -JA^BETA ) 

0, B,I 

1 RN=N 

A:^S=ALPHA-RrTA 
AB= ALPHA + BFT a 
ARN = ARfRN+ 1. 

A^!= ALPHA + PN 
8V=BETA+RM 
C3L=C?l=AMB 
C12=-R\i 

Cl l=C?l =^«N/( AN4-RM) 

" 1 R = 2 •«AN*Rf'I/( AN^-BM ) 

C?2=ABv2. 

C23=r-P\^ABNI 
C 32^C?2^-?. 

C33=C22+r23 

x:=-l.f2.«(r.ETAi-U) / (^3N+(RN“1.) *SQRT( AN*ABN/(3ETA+2.))) 
00 ^ K=1t23 

CALL JREC UP (N, ALPHA, BETA, XCfP, PI) 

P S= 

PSl=Pl 

xs = xr 

5 X= ll/I 1 .-XS^^'XS ) 

P3=((C11 + CL2*XS)^PS^”C13=#'PS1)<=RX 
P?3=( (C2 1+C2 2*XS)*P 3<>C2j^‘=>S ) *R X 
3=DS/P0' 

V=l •“U*P?n/PD 

H=-PNi^;^3/( l.+ S(vRT( (RM-l. )*( Dfv{«v-U ) ) ) 

XC=XC^-H 

IF(K.l T.3) on TO 2 

IF( H.L T. l.F*-24.nP.XC.F0.XS ) GO TO 3 

2 Cn\JTlMJF 

3 CIHTMJF 

CALL JRF C UP( N, ALPHA, RET A, <C,P, PI) 

X( 1 )^XC 

W ( I ) = CO'i* ( I . -XC*XC) / PI / PL 

AS=bS=D- 

3D 7 M=?,N 

P3 3= ( (C3 I +C 3 2*XS )^P 20+-C 3 3-PO )>»^RX 

0 = . S^'PZU/PO 

U=0-AS 

V= ( hJ+0-P3f)/P0/ i.-HS )/{ U*J) 

R M= 

R 1 = RM4-1 . 

MI = M-l 

XC = XC-R l/U/( 1 . ^-SOR T (RH’f' m * V-l . ) ) ) 

11 5 K=l,23 

CALL JRECLR(N,ALPHA,3ETA,XC,P,Pl) 

PS = P 

PSl=Pl 

XS=XC 

RX=l./( l.-XS«XS) 

P i3 = { ( C 1 1 +-C 1 2 =*X S ) *P S f C I 3 «P S I ) ’^R X 
P23= ( ( C2 H-C22*XS) *Pr)+C2 3^PS ) «RX 
AS=R5 = D . 

DO 4 I =l ,Ml 
0S= l./( XS-X( I ) ) 

AS= AS^-DS 
4 BS=RS+DS*DS 
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o=f>s/p:) 

U=l .-AS*0 

V=( l.-3^P20/P:l-RS«a*Q)/(U*U) 

H = -Rl=«‘3/U/( I . + SQRT( RM^(Rl<'V-l. ) ) ) 

XC=XC+H 

IF(K.LT.3) GO TO 5 

IF (H.LT , l.F-*24,nR,XC.EQ.XS ) GO TO 6 
CO^TIMUF 

6 CONTINJF 

CALL JR£CUR( M , ALPHA , 3E TA , XC f Pf PU 
X(M)=XC 

<< (M)==CON*( l.-XC^'XC) /^l/Pl 

7 CD'^TINOP 
RETURN 

8 Al = ALPHA+U 
A = A l + l. 

Bl=6ETA+l . 

C = A +B"~ 1 • 

S = OSORT ( AAB^'C) 

0 = 4.*CON*(A<*B+S) 

X=P,*R*R1/(B*C+S)-1. 

W = 0*B/(Bl«(S+B)*<'2) 

X ( 2 )=l.-2.*A«Al/( A^C^S) 

^ ( 2 ) = 0>«'A/ ( Al^{ 

RETURN 

9 IFIN.LT.URETURN 
Al = ALPHA4-1 * 

Bl = BFTA + U 
AB2=AlfB I 
A2= Al+1 . 

B2=R1+1. 

X=( BETA-ALPHA)/ AB2 
W=CON*( l.-X^X) 

W = W*2.*A02’^( 1. /Al + l./3U=«'<'2/( A2*B2>«'( I , / A 2<- L . /B 2 I * *2 ) 
RETURN 

END SUBROUTINE JGAUSS 



SU^RQllTlNf JRFCURCM, ALPHA, BETA, X, P, PI) 

CTMr'UTPS nth ANO (N“l)TH JACrv^I PGLYNnMiALS 

TOJWLE ALPHA,BETA,X ,P,Pl,P2, AM,B'^,CM,DH,G,C0,C^ 

AM= ALPHA 
3 H=BETA 
P2=0M=1 . 

CO = AM-i-BM 
r; = AM~BM 
CM = C04-l . 

= . 

P1=.5>!=(G + C?*X) 

G=3^C0 
DC) 1 H=2,M 
A'^ = A'14-1, 

BM=RM-H . 

f)M=DMi-l • 

CO = C? 

C? = C 2 ^? . 

CM=CM4-1 , 

P = ( CO*-l. )’?=(G+C0*C2*X) =«'‘Pl-AM*BM*C?4=P2) /(DM<'CM«CO) 

P2=Pl 
1 P1=P 
P1 = P2 
RETURN 

FNT SUBROUTtNF JRECUR 


SUBPQUrTNL JWTCON(C ,N, A,fl) 
nnUBLE C , A,B , AM,BM, ABN, DH, OGAMF 

COMPUTES THF constant t^AHT TF THF JACOB I-GAUSS QUADRATURE WEIGHT. 
AM= A^- 1 . 0 

A6H = AM4-BM-1. 

DM= L. 

C = .2 5*? .0*’»'ABM*DGAHF ( AHI^DGAMFIB^D^AM^'BM/IOGAMFIABM) =<^ABM) 

10 I M= 2 ,N 
AM= AM^-l . 

RM=BM4-1 . 

•A 8 M = A 8 H + l- 

1 C = C=«=AM/0M*B^/A6M 

C = l. /AM4-1./8M)^4=2 
RETURN 
END 


SUBROUTINE JWTCON 


SJRR'IUTINE LGAUSS(M ,XfWf ALPHA) 

COMPUTFS LAGJFPRF-GAUS^ Ai3S:iSSAS + WFIGHfS FOR ALPHA. GT. 
OOJRLF X( l)tW( I )t ALPHA 
^'UBLF XC ,PtPl fCON, Al, A2, ISQRT 
CALL L^TCON(CUN»N, ALPHA) 

I F( N-2) 9, 3, I 

1 XC= ( 1.4- ALPHA) / ( 1 .<-( N-1 . )/S0RT(2. + ALPHA ) ) 

RN = N 

AN=^ALPH A + Disj 

DO 2 K= l» 23 

GALL LRFCURIN, ALPHA, XC,P, PI ) 

PS = P 
PSI=P1 

xs = xr 

P:) = (RNi^PS-AM*PSl) / XS 

P2 9= ( ( XS-l A1 PHA)’J'PD-RN*PS)/XS 

0=PS/PO 

V=l.-0«=P20/PD 

/ ( l.+SORT ( (RN-1 . )«( ) ) ) 

XC = XC^-H 

IF(K.LT. n GO TO 2 

IF( H.lT . I.E-24.0R.XC.E9.XS) GO TO 3 

2 CONTINUE 

3 rONTINJF 

CALL L^ECURIN, ALPHA, X:,P, PI) 

X( I )=XC 

>^(1 ) = CON’<'XC/Pl/Pl 
A=R=0. 

DO 7 M=2fN 

P33 = ( (X$-2.> AL PHA)*P2D-(RVJ-l . )*P<)) /XS 

0=. 5^P20/PD 

U=0-A 

Q«0-P3O/PD/3.-B) /(J«U) 

«m=n-h 

Rl = RH4-l. 

M1=M-1 

XC = XC-R l/U/( 1.<'SQRT (RH*( Rl*V-l. ) ) ) 

DO 5 K=l,23 

CALL LRECLIP(N, ALPHA, X:,P, PI) 

PS=P 
P5l=Pl 
XS = XC 

PO= (RN^PS-AN’i^PSl )/XS 

P2T=( (XS-l.-AL PHA )*PD-RN’«'PS) /XS 

A=3=0. 

DO 4 1= l,Ml 

) = i ./(xs-xn ) ) 

A = A4-0 

4 r=B+D*0 
0=PS/PD 
U=l.-A*0 

V=( 1. -0* P2 0/ PD-B*Q* 0) / (U^J ) 

H=-Rl*Q/U/( I. + SQRKRM^IRH'V-I.) ) ) 

XC=XC^-H 

IFCK.LT.B) GO TO 5 

IF( H.LT . I.E-24.0R.XC.EO.XS) GO TO 6 

5 CONTINUE 

6 CONTINJE 

CALL LRFCURIN, ALPHAfXCf P,Pl) 

X( H)=XC 

/J(M) = C0N*XC/P1/PI 

7 CONTINUE 
RETURN 



8 A2=ALPHA+2. 

X(2)=A2^0S0RT( A2) 

W(2 ) = C3N*A2/X( 2) 

A1 = ALPHA + U 

X ( I ) = A2*A1/X (2 ) 

W(n = C3N*X(2)/Al 
RFTURN 

9 IF( N^LT^ 1) RETURN 
X = 1 .+-ALPHA 

^ = 2 .*C3N’f'( X^*l• )/X 
RETURN 

END SUBROUTINE LGAUSS 


SUBROUTINE L R F CIJP ( N t AL P HA . X , P , P 1 ) 

DOUBLE ALPHA,XfP»Dl ,P2tRM, AM 
COMPUTES NTH ANO (N-l)TH LAGUEKRE POLYNOMIALS 
RM=l. 

AM= ALPHA 
P2=l. 

PI = PM + AM-X 
DO 1 M=2,N 
KM= RM + l • 

AM=AM^l. 

P=( (R M4- AM-X ) 1-AM*P2 ) /RM 

P2 = P I 
I P1 = P 
P1 = P? 

RETURN 

END SUBROUTINE LREC’JR 


SUBROUTINE L WT CON ( C , M, A ) 

DOUBLE C, At AM,OM tDGAME 

CO'^PUTES THP CONSTANT PARI OF THE L AGUERRE-GAJSS QUADRATURE k^EIGHT. 
AM= A+U 0 
:=0GAM= ( AM) 

DM= !• 

DO I M=2tN 
C^C^AM/DM 
AM= AM^‘l . 

I OM=DMi-l. 

C=C/ ( AM*OM) 

RETURN 
ENO 


subroutine lwtcon 



mUBLE FUNCTION OGAMF(X) 

•VJJBl E R(7)»C(7) ,nLOG,DFXP,X,F,Z 

C1MPUT-S GAM'^A FUNCTION FJR :^EA^ ARdiMENT. DGAMF ANU X MUST BE TYPED DnUBLR IN 
CALLING PROGRAM. X MUST BE LESS THAN 144 AND NOT A NEGATIVE INTEGER. 

DATA F / Z.S0662 82 746 3U00 50241 57652 84810 0 / 

DATA B/l.n,l.n,5 3.J, L 9 5 . 0 1 2 2 999 . D , 2 99 4 5 2 3 . 1) t 1095 35 241 009. D/ 

DATA C/ 12. 0,30. 0,21 0.3, 371. Dt 22 73 7. 0,197331 42. D, 482 642 75462 . D / 

CHECK FOR X AN INTEG'^R. IF SO GO TO 3. 

N=I^INT(X) 

Z = X-N 

IFIZ.F'j.O.) GO TO 3 

CHECK FQO X A HALF INTEGER. IF SO GO TO 5. 

N=lDINT(X-.5) 

Z=X-.5-N 

IFIZ.EO.O.) GO TO 5 
N= 145-X 

I F ( N.LT. 1 ) STOP 3777 
Z = X+N 
DGAMF=0 
00 I 1=1,7 

1 3GAMF = 8( 8-1 ) /Ct 8-1 ) /(Z+OGAMF ) 

OGAMF = F*l)FXP (DGAMF-Z + ( Z-.5 ) +DLnG( Z) ) 

DO 2 I=1,M 

2 3GAMF=3GAMF/ (N-I+X) 

RFTURN 

3 3GAMF=1. 

Z=D. 

DO 4 I=2,N 

Z = Z^-l. 

4 DGAMF=DGAMF*Z 

RETURN 

5 nGAMF=i. 77245335090551602729816748330 
TF( N)6, 8,9 

6 N = -N 
Z = .5 

00 7 I=l,N 

Z=Z-1 . 

7 DGAMF=3GAMF/Z 

8 RFTURN 

9 Z=-.5 

00 10 1 = 1, N 
Z=Z+l. 

10 f)GAMF = OGAMF^Z 

RETURN 

END DOUBLE FUNCTION DGAMF 
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Illllllllllll 1 1 I 


nil lllll■lll 


1 1 I I nil II I I I II I I 


SnBQO'JTlNB Mr,AIJSS(N,X,W) 

COMPUTES ABSCISSAS + WEIGHTS FDR A NON-CLASSICAL GAUSSIAN OUADKATUPE 
CONSTANTS R AND C MUST RE FJRNISHFO BY USEP 
CnMMON/ROFN/R( 100)/C0FN/C( 100) 

U0J8LF 0,C 

DOUBLE Xm,W(l) 

nnjBLF XC » PN, PNl, PNP,^M2P» HN 

RN = N 

KNl=RN-l. 

P = -B 
RP = D. 

HN=C 

DO 2 M=?,N 
HN=MN*C(M) 

PP=rp-r«B(M)-C(M) 

2 R=R-B(M) 

XC = -(R4'5QRT(PNl=k(RM[*R4P-2.#KN=«'RP)))/^N 
DO 3 K= Lf 23 

CALL PNRECUr (r-UXC, PN,PN1 f PNP, PN2P ) 

X5= XC 

0=PN/pNP 

Y = l ,-Q^PNZP/P^lP 

3X = -RN* 0 / { 1. ^-SQRT { RNl«= <RN^V-1. ) ) ) 

XC = XC+-DX 

IFIK.LT.3) no TO 3 

I F( OX.LT . I .F -2A.nR . XC. EU. < S ) GO TO A 

3 CONTINUE 
A CONTINUE 

CALL PNRECU^ ( N , XC t P N , pn 1 , P NP , Pfsl2P ) 

xs=xc 

X(1 )=XC 

M ( 1 ) = HN/( PNl=f*PNP) 

AS=RS=0. 

DO S M=2tN 
RM=N-M 
Pl = RM4'l . 

Ml=M-l 

XC = XC-PNP/( .5*PN2P-AS-PNP) 

DO 6 K= I , 23 

C'aLL PN RECUR (N,XC,PN,PNl,PNP,PN2P) 

XS=XC 

AS=BS=0. 

DO S 1=1, Ml 
0=1 ./(XS-X( T ) ) 

AS = AS^-D 

5 BS=BS+D*D 
0=PN/PNP 
U=1 .-A5*0 

V= ( l.-Q*PN2P/PNP-BS#3*D ) /( U*U) 

DX=-Rl^Q/U/( l.+SQRT { RM* ( R1 « V-1 . ) ) ) 

XC=XC+DX 

(F(K.LT.3) GO TO 6 

IF(DX.LT.l.E-2A.nR.XC.EQ.<S) GO TO 7 

6 CONTINUE 

7 CONTINUE 

CALL PNRECUR (N , XC , P N , P N 1 , PNP , PN2P ) 

XS=XC 

X(M)=XC 

W(M)=HN/(PNl«PNP) 

3 CONTINUE 
RETURN 

END SUBROUTINE NGAUSS 
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■ I II III 


I 


OnJBLF FUNCTITN PNFUN(X,N) 

nnjRLF x,Pi,P5 

COM.MnN/HOFN/B( 100) /CJFNJ/C( 100) 
DOiJHLF B,C 
IF(N-l)4»3fl 
I Pl= 1. 

PMPUN=X-B( 1) 
m 7 M=2,N 
PS=PNFJN 

PNFUN= (X-rt (M) ) =^PNFJN-C ( M ) <^P 1 
7 . Pl = PS 
PETUPN 

3 PNFUN=X-B(l) 

^FTURN 

4 PNJFllN=l. 

RFTIJR*^ 

FNO OOU3LE FUNCTinN PNFUN 


SUSP our TNE PNRECJP ( N» X r PNf PNlt PNP»P\’2P ) 
CnMMON/BOFN/B( 100)/COFN/C( 100) 

03JBLE BtC 

OOURLE XtPN, PNl,PNPTP\2Pf XB, PNPlf PN2PI ,PS,PSP» PS2P 
PNP1=PN2P1=PN2P=0. 

PMl=PMP= 1. 

PN=X-B 

IF(N.lT*2) return 

00 1 M=2»N 

PS = PN 

PSP=PNP 

PS2P=PN2P 

XR=X“R( N ) 

PN2P=XB’«'PN2P-C{ M) *PN2Pl+2. ^‘PNP 
PN2P1=PS2P 

PNP = XB<'PNP-C(M)*PNP l^PN 
PNPl = P5 P 

PN=XB<^PN-C(M)*PN 1 
I PNl=PS 
RETURN 

FNO SUBROUTINE PNRECJR 


SUBROUTINE P NO ER ( Nt M, X t P3 , P $ ) 

OOURLF X ,PD( I ) »PS( 1 ) fXR ,S 

COMNON/BOFN/B( 100)/C0FN/C( 100) 

OGJBLF R,C 

* COMPUTES ( O/OX) **K*P(X )/(<-F4CT0R[AL) FO^ K=0 THRU M-1. RESULTS FOR 

* POLYNOMIAL OF DEGREE N ARE STORFO IN PD. RESULTS FOR POLYNOMIAL OF 

* DEGREE N-l ARE STORFO IN PS. M MUST BF AT LEAST 2. PD»PS MUST BE 

* DIMENSIONED M OR LARGER IN CALLING PROGRAM. N MUST BF AT LEAST 2. 

P0( l)=X-B(l) 

DO I K=2,M 
PS(K) =0. 

1 PD(K)=D. 

PD( 2)=PS( l)=l. 

00 2 L=2fN 
XB=X-B< L ) 

S=PD{ 1 ) 

P0( l)=XB*S-C(L I^PSI 1) 

PS< 1)=S 
DO 2 K=2fM 
S=PD(KI 

PD(K)=XB«S-C CL )*pS(K)4-PS(<-1) 

2 PSCK)=S 
RETURN 

END SUBROUTINE PNOER 
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SUBROUTINE CO^^PARE 

COMPARES 2N0 THRU 16TH RECURSION COEFFICIENTS WITH THOSE TABULATED IN HEF . 
COMMON/rtOFN/B(100)/COFN/C( 100) 

DOUBLE 3,C 

DOUBLE 8S(l<S)»CS(l6 ),ERR1R 
^EAD ID, BS(l),CSm 
IF( EDF,5)4,l 

1 BMAX=CMAX=0. 

MB=MC=0 

no 3 M=2,16 

READ ID, RS(M) ,CS(M) 

TFST=DABS(RS(M)-B{M)) 

IF( TEST.LF.RMAX) GO TD 2 

RMAX=TFST 

MB=M 

2 T^ST= )ABS( CS ( M)-CI M ) ) 

IF( TFST. LF.CMAX ) 30 TO 3 
CMAX=TEST 

MC = M 

3 CONTINUE 
EPR0R=3(MB)'PS(MR) 

PRINT 11, MB, B( MB) , MB, 3S( MB ) ,ERROP 
ERR^R = C(MCHCS(MC) 

PRINT 12, MC , C ( MC) , MC, CS( MC ) .ERROR 

4 RETURN 

10 RORMATI ?D40. 30) 

11 F0RMAT(//4H. -B(,I2,3H) =,036.28/ 

f 4H BS( ,12, 3H» =,036.28/ 

f 9H ERROR =,Q36.28/) 

12 FORMAT(//'+H C(,,I2,3H) =,036.28/ 

¥ 4H CS ( , 12, 3H) =,036.28/ 

+ QH ERROR =,036.28/) 

ENO SURROUTI NE COMPARE 



Sample -Problem Program Input 


For the case ALPHA = 1 the program input consists of six data cards containing 
the following parameters with formats included in parentheses: 

IFLAG(9I4) 

ALPHA, UA(2D40. 30) 

LMAX(I4) 

NC(I4) 

LX,NMAX(2I4) 

N1,N2,N3(3I4) 


The input to the sample case is listed below as card images. 


000000000111111111122 22 2 22222 35333333334444 Mi 444^^55 
12345C789012345C7CQ01254SC7C901234S67C9012345C7C901 


5 35 35 55 3C6C6C666C6 77777777778 
23430719012343678901234567890 


0 

1.5 

C 

1400 

3 96 

8 96 


I 


.613705638880109381165335757 


UOOC 000001111111111 22222222 22 33 333333334444444444 3 33353353360666666 66 7 777777 7 778 
12 34 36 789012 34 36 78901 234 36 789C1 234 36 789C1 23436 78901 2 34 30 789C1 2 34 36 7890 12 34 56 7890 


Sample -Problem Program Output 

The sample program output consists of a line printer listing of the correction coef- 
ficients (see eq. (37)), the recursion coefficients bn and Cn (including hQ 

stored in cj) punched on cards, and the requested quadrature abscissas and weights 
punched on cards. Tables I, n, and HI were copied from the sample -case output files. 
The quantities output in tables I, II, and III are presented using a 29 -significant -figure 
format because this is the double precision word length of the computer used. However, 
the numbers in tables II and III are only accurate to about 18 significant figures for the 
sample case. 


TABLE I. - CORRECTION COEFFICIENTS AND THEIR LOGARITHMS 


CO^RPCTION COEFFICIENTS FOR ALPHA = 1.500000000000000)0000000000000+00 

AND NC = 1400 


CIO) 
C( 1 ) 
C(2) 
C(3) 
C(4) 
C( 5) 


3.21789035155754345691211094010-07 
-6. 52323042044022828056934274850-14 
4.38187106178325216340494404670-20 
-5. 971 I 52693826043220 5944689451 0-26 
3.886991445 17711611173387556700-28 
2. 36406220651858388562855804450-28 


LOGIOI ARS(CtO) ) ) 
L0G10(ABS(Cmn 
LOGIOI ARS(C(2))) 
LOGIOI ABS(CI3 ) ) ) 
LOGIOI ABS(C(4) ) ) 
LOGlOl ABS(C(5) ) ) 


-6.49 

-13.19 

-19.36 

-25.22 

-27.41 

-27.63 


37 


TABLE n.- RECURSION COEFFICIENTS 

‘^eCURSION C3EFFIC1PNTS FUR r HF WEIGHT FUNCTION L N ( I . /X ) <■ ( 1 . - X » ** ( ALP H!\-l . ) OVER 
THE INTERVAL (0.,l.). 

ALPHA = 1.50000000000000000000000000000«-00 

LX = 3 

NC = lAOO 
LNAX = <56 

N BIN) CIN) 


1 2. 12554 52108 71229 350 309926 93 040-01 3.53 58 153 7031 1 8403 18881 349490 7D-01 

2 4.333431 11361 125328989241440530-01 3 . 8026 1 22 53094 772034470 1 4039220-02 

3 4.70516450471843969719383403740-01 5.3551 51 1 959857 1 1 3 5 7145 71 5472 10-02 

4 4. 83365 62 32 7826 09660 52 58 8 7 79460-01 5. 79 1 0 73 3 3 8 7526 0009 7907 7283 1 3 OD- 02 

5 4.89313473482 8178 32642 33 53 1 7320-01 5. 9718 866 6573 568197 709400309960-02 

6 4.92552923213164195250293412570-01 6. 0637 353642 3 1 0478 7034299822890- 02 

7 4.9451218C89 5151 1 128 6052 4374700-01 6 . 1 1 66488 1 09233 09092 073 7700 38 20-02 

8 4. 957 8 75 3649 29529 262 704 8 72 70740-01 6. 1 498 7054883868407 2592 1 590 22 50-02 

9 4.96664241 710 3399 976 5 7949 491290-01 6. 1 720 78 34 14353 92 1 6960785 1 52320-02 

10 4.9 72 9 2 83 5 392 2718 023 0595 77 02930-01 6. 1 876494 61 836 5 94283360 704635 50-02 


1 I 4. 97758909664209024456152654620-0 

12 4.981 14068445184391387401460990-0 

13 4.983909 34497 7413 755 51836199 740-0 

14 4.98610954917032304031371800140-0 

15 4.98783701 575961047380358974770-0 

16 4. 93934357030356653499512493920-0 

17 4.99055211 279701362061908089130-0 

1 8 4.991565949136437421302511 073 OD-0 

19 4.99242479029029405367481883280-0 

20 4. 99315371 153180062235156553400-0 

21 4.993 790*82 9061 199 672 6949 4346750-0 

22 4.99433914681007707766099942800-0 

23 4. 994 81785 631 7769 72 75 32301 39980-0 

24 4.99523826613743466043223330600-0 

25 4. 99560947736645927736374203790-0 

26 4.99593383139489128628066671850-0 

27 4.9962325314382665001341 9252320-0 

2 8 4. 996 4 9542 31 5 5835 064 7949 7563480-0 

29 4.99673170892422282914454411050-0 

30 4. 996 94486 3 1078 12 882 8 1 50 Cl 37030-0 

31 4.99713781073120166880236388000-0 

32 4.99731302853151063638338399980-0 

33 4. 99747262496616057195569276830-0 

34 4.99761840404039668813534631 790-0 

35 4.99775191 6589623 752 0456 9478590-0 

36 4.99787450175786238092473672040-0 

37 4. 99798732075761020383535866260-0 

38 4.99809138451027784531716901640-0 

39 4.9981875764030341395381 7876740-0 

40 4.9982766711 2 39B04364277 7324730-0 

41 4.99835935032954529678345813110-0 

42 4.99843 62157 387542 322236 9672 000-0 

43 4.99850780012631613866156479140-0 


1 6.19898539395190034161136539360-02 
1 6. 20749254381854953955989954250-02 
1 6.2140387899955068 1140835105350-02 
1 6.21918315988199200833937987310-02 
1 6.22 32 9 88 73 86 733645 0340385415 70-02 
1 6. 22664279042924204879611806830-02 
1 6. 2293 9636093231910743548567240-02 
1 6,23169070717321311190149133450-02 
1 6. 23362248436493690975062 1905 10-02 
1 6.23 52 64> 046586 71532 5929462 78 10-02 
1 6.23667112614443300100487115560-02 
1 6.23783595350580233065834478300-02 
1 5.2389421 1003 974524 73 13001 309 70-02 
1 6. 23986605646331609387073107110-02 
I 6.24067896109647849645463221050-02 
1 6.24139792074903113714400962180-02 
I 6. 242036 865293 7 784509 594453 6490-02 
1 6.2426 0 72 3621716 7199043 352 752 70-02 
1 6. 2431 1850145113864870323003550-02 
1 6.24357855011423320233526145560-02 
1 6.24399399813694339653053648950-02 
1 6.2443 704 2703 99768880 32 616015 70-02 
1 6.24471257206502212471044332060-02 
1 6. 2450 244 71571795525135 0813 78 70-02 
1 6.24530958655302318511221737260-02 
1 6. 24 55 7089690694765 8917054131 00-02 
1 6.24581097949242000814283067660-02 
1 6.24603207180178537868659677810-02 
1 6.2462 36124201887217439 32241760-02 
I 6. 24642484302977865497273853390-02 
1 6. 24659972632783 502967490151050-02 
1 6.24676209362053891795862843010-02 
1 6.24691311084170642460010556870-02 
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TABLE n. - Concluded 


4^ 4. 9985 74576 591 2 64944 06 95 8 882 80ri-0 1 

45 4.998 6 3696640 3 561 596 I 6 1 8 44 31 460-0 I 

46 4.99869534567267701 133944598560-01 

47 4.993 75 00*^1036054133 3 560 51475 70-01 

48 4.993 801 3 845 2871 2B906 I 52 05 821 10-01 

49 4.99384961 776600963648892192250-01 

50 4 .99889499554809374596300656500-01 

51 4.99893773897567624367600885740-01 

5? 4.9989700481514004933 1220116760-01 

53 4.9990161045286346571 1824522490-01 

54 4.99905207295932005081223994480-01 

55 4.99908610 348415 48840 55 3 7047 380-01 

5 6 4.9991 183*32^015138621996 8040150-01 

57 4.9991 4888614581 896774510664810-01 

58 4.9991 77877501361 401 7160 834529D- 01 

59 4.99920541 167365020553487252880-01 

60 9992 3 15R 47 3 7035 98744794 591090-0 I 

61 4.99925648497501287925914797550-01 

6? 4.99928019362589493137069895860-01 

63 4.999302735547^1882064370968270-01 

64 4. 999324329808681 27606825394190-01 

65 4.999344360? 1922658837572036160-01 

66 4. 99936452580251561347799998640-01 

67 '+.9993 83 29 12 2039 8 794 2 450 1 772 380-01 

66 4 .9994012371 543 03 84 8 8 05 2 5 72? 490 -01 

69 4. 99941 841 06481 I 06 79 240^6209890-0 1 

70 4.99943485*^417052 463 3 719 1 161270-01 

71 '+.9994 50612126 43 6 75 48385 72 554 70-01 

72 4.99946571 864351 03 360829 33 74490-0 1 

73 '+. 999480210265 38 5 329 3 53 6 85 96 820-01 

74 ^.99949411 992559226798154983690-01 

75 4.99950747838 152053637458095210-01 

76 4.909520314384741 13 77059 29082 10-01 

7 7 999 5 32 6 548 3 597548764708887350-01 

78 4. 99 95 445 24926 27 2 08 3 19592 068o50-01 

7 9 '+.999 5 55 9482 6 5 776 3 92 88 42 9 8 643 50-01 

8C 4.999566947001324 70092 684440180-01 

Bl 4.999 5 7 754192 3956 96 3 5 99 3 5 8 2 5240-0 1 

82 4.99958775256732448004740152500-01 

33 4. 99959759729786313241159366160-01 

84 4.99960709339751075832903230490-01 

85 4.999616257! 3966403902922554720-01 

86 4. 99962 51 03 85 9000 745 5436 2 3 52 32D-0 I 

87 4.99963364801572559926536824160-01 

88 4. 99964190325474394656951332410-01 

89 4.99964988246021553382092977120-01 

90 4.99965759780589612799387534140-01 

91 4. 999665060801634641633551 84120-01 

92 4. 99967228233635770561029704310-01 

93 4.999679272717841 732 0090 73 2 3420-0 1 

94 4.9996860417 09543 973 99 35 9 942310-01 

95 4.99969259856473858199387898180-01 

96 4.99969895205813061150636813870-01 


6.2470 5 33 11294117482 69862 75 7980-02 
6. 2471851 13347103965568198242 10-02 
6.24730783543958852 175200634160-02 
6. 24742270884738059766951400070-0? 
6. 24753038858749330381649089010-02 
6. 24763146275376279124536409010-02 
6. 24772646060326746277605221760-02 
6.24781585938900331779078826080-02 
6.247900090458701376804102291 10-02 
6.24797954441031532352720649720-02 
6.24805457563519478107763437380-02 
6. 24812550626672289320829622290-02 
6.2481 9262964439609296333703050-02 
6. 24825621334529594127352316940-02 
6. 24 83 16 5018 5212175107898795780-02 
6. 2483 73 7 1 8903999403928 1 32 52 8 80- 02 
6. 24842805957282114306067522310-02 
6.248479742101 13395968715584070-02 
6.24852890953236336298543926200-02 
6.24 85 75 73115 9759085730 35 313710-02 
5.248620353816736831963718441 30-02 
6. 24866291302314909288464401290-02 
6. 2'+ 87 03 5 340 393 5 3087 106 496124 30-02 
6. 248742332737676830612 129941 80-02 
6.248779^1647895082059720159340-02 
6.248814884830119637160041471 80-02 
6. 2483 433 30 23 753 085 9019 702 030*90-02 
6.24888133862928 143576894694650-02 
6. 2489 12 48995 8953545940656606 10-02 
6. 24 894? 3 5 8697 1 6 75 01 073 70 1 5 00 ID- 02 
6.2489 710142 765 970 7 0432 37 70 8440-0? 
6.2489 985 214954 147 57 203 49 56 03 80-02 
6. 24902494088354594681023409620-02 
6. 24905032903559871863729535300-02 
6.24907473891388967376534323460-02 
6.24909822012459623263821981660-02 
6.2491 208 19169724720373879411 10-02 
6.2491425 796 7 72 8 455 8754 9509 79 00-02 
6. 24916354261179636772183323550-02 
6.24918 3 74646 703 09 08264 96 27 92 00-02 
6.24920322744267242708044909190-02 
6^ 24 92 22 0196 0642 05 3 0 8459926 63 50-02 
6.2492401 550428861611 1014463310-02 
6. 24925766399049691888638162380-02 

6.2492 745 749675 02620658 72 950400-02 
6.249290914883061581791 766585 70-02 

6. 2493 06 7091692 9000960204 83 4730-02 
6.24932198183006956244455571020-02 
6. 249336755582330302615904034 70-02 
6.2^935105191545681577887448790-02 
6.24936489117440320761081924670-02 
6.2493 782 926 32 04716 3003 32 54 38 30-02 
6. 249 39 I 2 7455626 35 14996 33 02 45 80- 02 


39 


TABLE in. - GAUSSIAN QUADRATURE ABSCISSAS AND WEIGHTS 


ABSCISSAS AMO WEIGHTS FOR THE WEIGHT FUNCTION LN ( I ./X ) *( I . -X I ♦* ( ALPHA-l . ) OVER 
THE INTERVAL FOR ALPHA = I . S 000000000 00 00 000 000000000000+0 0 AND 

FOR N = a THRU N = 96 IN INCREMENTS OF 3 


N = 8 


M XI M ) 

1 1 .2680038934'+52 11612784557 81020-02 

2 7. 5901 1957704154566324203934700-02 

3 1 .8 85 62 3166 3 3 5672 36 946192 6 04010-01 

4 3 .3879 53 07 3? 3 489 3565 6 743 04 765 10- 01 

5 5. 089398949^14471 37368963391340-01 

6 6 . 789 5356659561 1 7 497R8 1641 12000-01 

7 3.23482351400709740743624243520-01 

8 9. 39530648627747992141395893120-01 


W(M) 

1. 5718 2031727 7.86014 9814 082867 20-01 
2.21 749345279989400364783489760-01 
2. 01280722783348405870128523340-01 
1. 43764455232332422352546552430-01 
8. 21 38641 5474023900201612512370-02 
3. 5772 063 912403 149 83 941496948 00-02 
1.03616949254898089275360331540-02 
1. 33258162243243953215534348280-03 


N = IS 


M X( M) 

I 3. 793 1 1544977950313470036936350-03 

? 2. 24023555072729901234335013930-02 

3 5. 6 702 3769 34 2 2^247 373 9846 8 131 5 10-02 

4 I . 05777629534980422151556985030-01 

5 1 .681052690472364920691387801 lO-O I 

6 2.41664686129631752720069327140-01 

7 3.24023911 53231 3867193592158140-01 

8 4.1243 063815621911 363 I 76 7 716480-01 

9 5. 039 10 325 I 7919 73 708 6 0103 316420-01 

10 5.9537342372 1555769579322961340-01 

1 1 6. 8371342052 80417 534412 94418800-01 

1 2 7.659 44 72 0 39 95 840 76 39 2 746 3 72650-01 

13 8.3975450071 1999312812223899160-01 

14 9. 01150403 84073492 835 336 033 9460-01 

15 9.495 2 2659 544 2259 52485137395210-01 

16 9 .82725382307036600369138116610-01 


W(M) 

5.931672016284623979031 12082320-02 
9. 9709 20 7 346 7904 1 745 785 822 523 7D-02 
1. 16778984846344003859961 178350-01 
1. 18959209377791045572076043500-01 

1. nil I 12 9446602 062 9058 356 64400-01 
9.69542 92 1 9 1 ? 77 1 92 9 5 5 70 1 1 5A 37 80-02 
7.955631492856140091 18001257680-02 
6. 1408 7544796 4900495298 5110 6910-02 
4.4405947067203581 1607907599070-02 

2. 9824686 3253304191868213078120-02 
1. 83402975942387187857398152300-02 
1.00894783406156620109111009410-02 

4. 7769154 760752 03 0 58545 779063 00-03 
1.81453297151413828921792009990-03 
4. 77347516430636483190206416890-04 

5. 7718 95911430450588935 05 848660-0 5 


N = 2'+ 




X( M) 


W(M) 


1 I .810912904143091 231 1 3383479570-03 

2 I . 0592036167781^209314266969790-02 

i 2 .68401 548 62 65065 986 2024 2 928 070 -02 

4 5. 0 38 7 3433 3 2 72612 79^8694944465 0-02 

5 8.001 330R049570469234312 7121050-02 

6 I . 17972443485372040061422602120-01 

7 I . 610074542710672626877959641 ID-Ol 

8 2.0936095615 1665 2 8 3 72647 1 3 61 5 00- 0 I 

Q 2. 62287358259784831081165172060-01 

10 3. 189 65 539 79 39068 ns 8 109 3 393930-0 1 

11 3.78512419848599980662731273860-01 

12 4,39997341515597053873160998660-01 

13 5.02457111614009279840762533680-01 

14 5.6491 1495405493876152078944490-01 

1 5 6.26378942183970167046761633830-01 

16 6. 85892305049594630644135310350-01 

1 7 7. 425 143131 5965 70 292 2666 2 953910-01 

18 7.95352555724756498070689099480-01 

19 8. 435 7 3 743 2792 504 5719 893 4803130-01 

2 0 8. 86 4 I 70 2 3 33 15342 6 32 8 560 7061960-01 


3. 1843 708 3666 00495 72 09 7 344768 60-02 
5.66662645507311420743097483530-02 
7. 1 1643473135756774879234724330-02 
7.89114854309674826816906828300-02 
8. 165962106380689923172540141 10-02 
8.05895464074182954017272900210-02 
7.66350640383285763687901859950-02 
7. 05 92 173162104 77912 6951 7106760-02 
6. 31 579795134616489215564187300-02 
5. 49424159642883362619267016220-02 
4. 6469305 1848758727935872980350-02 
3.81738669892935560086592461120-02 
3. 04 00 372 7 83 2 39 707 596441095 9360-02 
2.3401 7688594272954831659872770-02 
1. 7342 1809292098633461867341460-02 
1.230261 94665075371233030302600-02 
8.28980435435452276366192167400-03 
5. 2477 3 8019 34963 0087142 5885 86 00-03 
3.07 144572369652850327940768870-03 
1.6221 5525469218100478794876970-03 


40 


TABLE m.- Continued 


21 9.2 32 06146 31 12563 73 L 4653 8 I 37 093-01 

22 9. 53360314920977331631537399150-01 

23 9.76403683208362949142966314980-01 

24 9.91974631905422105086125424220-01 


7. 43099062066056491800609688310-04 
2. 75086071573599435267671393870-04 
7. 10046641311438274099310881070-05 
8.47880748320414586641824702940-06 


N = 32 


M X(H) 

1 1.0602 I 77045 66 333 003 55 32 3321550-03 

2 6.16029698036034040956995858200-03 

3 I .55988115738796185954024259410-02 

4 2.93325089132415292025344713290-02 

5 4. 72507330174069644589854453080-02 

6 6.92376947832017945751869129000-02 

7 9.50770219330400178532440565650-02 

8 1.24554337053324858041780116520-01 

9 l.57^f 12 5098807232 1080195 5621240-01 

10 1 .933 594254300882 5554 885 7994050-01 

1 1 2. 3? 0 75 88 890 1 82 33 5 17 803^» 7643680-0 1 

12 2.732 162085884428 14313000574660-01 

13 3. 1641 1963 4642364 792 1569 76392 30-01 

14 3. 6 127 5 444 97267 18 023 0531 8962340-01 

15 4.07403252995240685001348661540-01 

16 4.54380020226027798665436422700-01 

17 5.01782235478305325862074837240-01 

18 5.49132134020955480084679768990-01 

19 5.961516214678222515372591 77860-01 

20 6.42266196312989729290355420290-01 

21 6. 871 08839796994380744627133020-01 

22 7.30273827832829580308038647790-01 

23 7.71370443172413488142722374590-01 

24 8. 100 2 65 71 82654000851081 2618260-01 

25 8.4589205695 1145793104107146750-01 

26 8.78641950991377584770409708070-01 

27 9. 07979462911955463172166505240-01 

28 9.33638679656860882907570196190-01 

29 9. 553870058991544586l288607334n-bl 

30 9.73027319759588787836999676800-01 

31 9.86399904276508334301735 199310-01 

32 9. 95384869928502686597555318070-01 


W(M) 

2.01262612879950589456595715700-02 
3. 68892348994588273617289912180-02 
4.78795592302463302710527168010-02 
5. 51150060869643828016480848230-02 
5.95231897679062230526877882730-02 
6. 1678375706589944825667661 0380-02 
6.20015209180505578989401 101820-02 
6.0834457700744099221916264 7350-02 
5.8471972 646425 379630330091 2 770-02 
5. 51766217003230079354868808850-02 
5. 11855395223808352411719617560-02 
4. 67132748093130712399887485420-02 
4. 19526220333780193158751046360-02 
3.70745023937701336316892695080-02 
3.22274958247031351700765520640-02 
2. 75373841382190620259020786310-02 
2. 31069247283689565769055461010-02 
1. 901 59857565703652931239740190-02 
1. 53221 138064773344515546643980-02 
1. 2061 5617565510133200569098780-02 
9.25077136947194210511649334030-03 
6. 88828353093727530114395075870-03 
4.95 701439897542857824543105440-03 
3.42635197821472567183165902260-03 
2.25746678018217950234^80525260-03 
1.40127796842948741553043350640-03 
8.06462522650979099034996617530-04 
4. 19925333943971869059373266430-04 
1.90135105858514674652871760780-04 
6.97392977658976195861045561240-05 
1. 78770836926051313828293405830-05 
2. 12 5146625281.708743029766446 80-06 


N = 40 


M X( M) 

1 6.96^86704355641477354209026600-04 

2 4. 027 3 424 ’36 9055 18 774639 1506320-03 

3 1.01880346922303843? 14295708880-02 

4 1 .91663239915549725269294660510-02 

5 3. 092 2 65652392927 8795903 76 0713 0-02 

6 4.53971307500975820719808536310-02 

7 6.2 51 I 727979 59487 965 29836 203300-02 

8 8.217157250948283 16399630288280-02 

9 1.042 6 5903 850 106 190 1 7 3974926760-01 

10 1.286689686057137 17374220884300-01 

11 1. 5524091130103578636001 8496520-01 

12 1. 83828697927965516990608869830-01 

13 2. 14267086783070099712363670790-01 

14 2.46379651972476971229018605960-01 

15 2 .79979860090180997057270046710-01 

16 3. 14872197766425350007615492220-01 

17 3.50853346039695355834987758740-01 

18 3.8771 3396356659487614701840730-01 


WIM) 

1.3983238548631 1468766816300590-02 
2.61127532468732749658^48569820-02 
3.45619835415016044044727480830-02 
4. 06 36 640122 7069 899 33 701855 63 50-02 
4.492167786413799285376925671 lD-02 
4. 77636705430447373505986735750-02 
4.94037995765613018089545942750-02 
5.00285533013433860810649345780-02 
4.97927670435842467038329965060-02 
4. 883121177104953135514969146 10-02 
4. 72647572640378870768577084460-02 
4. 52037173986943852652122083590-02 
4.274963844132153973768 19122360-02 
3.99961949247254467341036 790300-02 
3. 70295685997479250003195488530-02 
3. 39285351022642712143592 152540-02 
3. 07643991356671185764433292500-02 
2.76008696621607529802241831900-02 


41 


TABLE m.- Continued 


19 ^.25237 l022l0104255D670155425in-Ol 

20 4.6320515985832935963822634229D-01 

21 5. 01 3 955 11 15 987 30071 485066 61 43 D-01 

22 5.39584661263986274628386142470-01 

23 5. 7754900369566701322902 1781840-01 

24 6.153661452491656 14725350037350-01 

25 6.5191622304621 3336499443718320-01 

26 6. 87383206118307369298767360350-01 

27 7.22756173931819395626180652350-01 

28 7. 56330564387183884246971561040-01 

29 7.88409383986491111031824154090-01 

30 8.1880437^075642654623905015120-01 

31 8.47337119326918630765455926140-01 

32 8.73340112936336704707333372190-01 

33 8.98157737373346749205569632960-01 

34 9 .201471 899577299034^8265242190-01 

35 9.396793271130006 19784775399860-01 

36 9. 56639430075429560448280986130-01 

37 9.70927888994143317783327241650-01 

38 9.82460811020855679403453797810-01 

39 9.911 70600278588027520739710190-01 

40 9.97006980881764682504174674340-01 


2. 44939360197186568572367611000-02 
2. 1491 78592 306846 189583576325 3D-02 
1.8634 7926126379550996096239740-02 
1. 59555885602281297772674684380-02 
1.34792357324005264272030684250-02 
1. 122349663381 17693573163222780-02 
9. 19920571306205243284500737090-03 
7.41073691609915296108761662300-03 
5. 85656003551898199100961551390-03 
4. 52 98 75 91719889809933510 53 92 80-03 
3.41931843035201296476477310450-03 
2. 5 097 098409917540 9 0839 77175810-03 
1.78285357482573716434208577770-03 
1.21835029497140749732498156270-03 
7. 9441 92376693592500378403521 30-04 
4. 88709431349962368777555369160-04 
2. 79084255979776285445460160450-04 
1.44 36 376 336469563 3 3 80087306180-04 
6.50100048740830941973533932180-05 
2.37417774420176293462301808280-05 
6.06662499766365184387069009860-06 
7. 1959 2567338439814565954 797220-07 


M = 48 


M X( M) 

1 4. 92770594991399739269654953150-04 

2 2.83378293739246880472285044020-03 

3 7. 17475677025330320560320373430-03 

4 1 . 34979619344721218038785757290-02 

5 2.1 7905470217157458844909705030-02 

6 3. 20240104441081382563081355930-02 

7 4.41605257858448384043119236050-02 

8 5.3153609^905323556319358465670-02 

9 7. 39^857513290371 70672553066420-02 

10 9. 14829133690217506815166782890-02 

11 1. 10686643 144069890535105204650-01 

12 1 .31432664583663365895426714130-01 

13 1. 53787120653000274694251825100-01 

14 I . 7 750977454 33839503 1931 229454, )-0 I 

15 7.02554405 182171454954425291240-01 

16 2. ? 83 19 221 8029063 95690930846630-01 

17 2 .5619 72 97 5194673 103 9266 2100240-01 

IP, 2. 845 7 7 021216 7921939 74242603 500-01 

19 3.1384256664861912859521 3358480-01 

20 3.438 74 37 733 6 5792 876 657109422 50-01 

21 3.74549665649695256244337753140-01 

22 4.05742924282 120642319355745140-01 

23 4. 3 7 32 6448 22 39338 5 1496 I 5 C572 030-01 

24 4.691 70865225696625171446698520-01 

25 5.011 45672683891701011046844780-01 

2 6 5. 331 1 9773 02 32814 44801885 345730-01 

27 5.64962046010548463228189630120-01 

28 5. 96541834012107928203804829820-01 

29 6.27729540995979319378514956530-01 

30 6.58397140073180362773597207340-01 

31 6. 884 18709225340338515976098490-01 

32 7.1767095234710067622966 1368150-01 

33 7. 46033709443898677674986536340-01 

34 r. 7339045371 7550306753468929370-01 

35 7.99628773603346834812964627590-01 

36 8. 24640837742905199920052427510-01 

37 8.48323840986615465415038200670-01 

38 8. 7058042959328595682B574672 530-01 


W(M) 

L. 03 3604447 15 783 85 8 5994 72 6038 80- 02 
1.95571906326880076124006107400-02 
2. 62283216289513050650682362770-02 
3. 12658012083979295729566854440-02 
3. 50756516281379729440776121070-02 
3. 7393 0310639244518053907041710-02 
3.98 76 758174095 7592 73 72 8002 3 700-02 
4- a 45147562252 5929 3384 20 73 5840-02 
4. 1793331 3312845136864617393620-02 
4. 19017281424005734470747686270-02 
4. 15409532210118121980949893420-02 
4.077465^6966591249167076687400-02 
3.96612037963313564879148855960-02 
3. 82547237034289617986624036650-02 
3. 6605 68 2 1780 9640476 573 05 777710-02 
3.47612371152886031534425070280-02 
3. 27654098155756902697543760450-02 
3.06591548916008151723810836220-02 
2.84803672289450787231298630030-02 
2. 62638544946391398662015786490-02 
2.40412952274593139754516807370-02 
2. 18 41 19680474915716712 330612 70-02 
1.96888635465144958670816617830-02 
1.76063822824081935192955223530-02 
1. 56126304989952523794562065510-02 
1. 37233104703450923829504733040-02 
1. 19510114059299810316441345820-02 
1.03053005298791965449455964240-02 
8. 79284307196195928559650474900-03 
7. 41755036179783328611304030800-03 
6.18075454626012565694149110400-03 
5.08140787761179029184194380230-03 
4.11630403446515009657922166830-03 
3.280318531091615198321 68337970-03 
2.56666493746526840062557172060-03 
1.96716336892097017153399758140-03 
1. 4725 175076 7232897 732641 3263 8D-03 
1.07259628644476892506975439360-03 


42 



TABLE ni.~ Continued 


39 8.913 19 10388309 L793357896862000-‘01 

40 9 .L0454596T00063Z46701 334640250-01 

41 9.279082261 5736894034271 1 718980-01 

42 9.43 60 3 221 28 74818 551*06845241320-01 

43 9. 5 749 00 19445 3608 2956 731 8894840-01 

44 9.69496534770302875989375377790-01 

45 9. 79578398312692854077857043710-01 

46 9.876941755303541 80299485762930-01 

47 9.93810596596319793715835327610-01 

48 9.97903129919343811855572936890-01 


7. 56 71 6295 133718894644045909210-04 
5.13920962691010157197334421110-04 
3. 33252616201949368358669295740-04 
2. 04013624896137605334722456610-04 
1. 16012995849160995350365916150-04 
5. 9794997209176431 1 7 149 322698 30-05 
2. 68466403522763582 180662220790-05 
9. 78 11 5047174730167637104389520-06 
2. 49489145620253264051501993010-06 
2.95 50166 16695326775 32538970010-07 


= 56 


M X( H) 

1 3.67181630506395381274704869210-04 

2 2. 10897478341050126212524295990-03 

3 5. 3258996 10220295 1125322 421 5360-03 

4 1.00184365022069254056159616750-02 

5 1.61 777 7 1<^ 16 1745 3 3443 4502 5 05 720-02 

6 2.3738896309951 8635687880366840-02 

7 3.2 83 I 48 142 I 6.?? 45 772 5289 7 5 726 50 -02 

8 4. 32802927818650100650391357930-02 

9 5.51054460644182472744028359090-02 

10 6. 827 >605 I 59 393 33 73 808 81 2 074020-02 

1 1 8.274814587653031 47594270367470-02 

12 9.84743227818398245349193720570-02 

13 I . 154 19426 369436 54183931 36 00810-01 

14 1.33 52 794 7486 88 59 950 0072 3 845260-01 

15 1. 52745747140433323613130853430-01 

16 1.7301 523351 1776662151041 772770-01 

17 I ,94275546768606000275866038000-01 

18 2. 164627S1 85 76802 5 32 9 75 8 7 72 7550-0 I 

19 2 .395 I 0 039 21? 41435 3532 6349964450-01 

20 2. 63347933067563533460405715950-01 

21 2.87904507 1785141 722 8769 7623670-0 1 

2 2 3.1 31 0 5607 I 88 2396 284 52 18 085 12 70-01 

23 3.38875079394533958251290762360-01 

24 3.651 35005153277763069637214970-01 

2 5 3. 918 05940 34224714892 082 6 72 91 00-01 

26 4.1830715921 1 3242 6C0 5 699 7 865 090-0 1 

2 7 4. 4 605 6 90192 3 2 590 863 7540 94 6 3370-01 

28 4. 73472625070925690592945859910-01 

29 5.0097125443713864901262 1930440-01 

30 5.23469439251460244192544341070-01 

31 5 .55833807 IB 86 2 00 561 49 68 3 55 I 850-01 

32 5. 8313 12 19347135865682 196194910-01 

33 6 . 101290244437453875621723281 80-01 

34 6.367953 1 1459223788265116 262730-01 

35 5.53049159973894901032280429360-01 

36 6 .888108874366566 14846142977740-01 

37 7. 1400229261 9331422512863585160-01 

38 7.38546894518792645809533242550-01 

39 7.62370165982321685853420507970-01 

40 7.85399761347073527901450670300-01 

41 8. 0756573740225905341721 1969870-01 

42 3.28800767002649574081437300490-01 

43 8.490403446842332 15797844461920-01 

44 8.68222983657302071300173512960-01 

45 8 .86290403578996514925143837830-01 

46 9.031 87708536619871982289268320-01 

47 9. 18863554705494378642692975020-01 

48 9. 33270307182347423171630576660-01 

49 9.45364185541297166262863430030-01 

50 9.53105397725705380659081637260-01 

51 9.58458262009274263544153138620-01 


W(M) 

7. 98 18624930081974797 1809941420-03 

1. 52 53024 604093 8 71302974932 3930-02 
2.06530657469156827882021949930-02 

2. 4861221 7602669598608910553930-02 
?. 81 769890941233127694497793590-02 
3.0771 8057860426072388780782640-02 
3. 27592208022380054867222854500-02 
3.42220824416665284010806462980-02 
3. 5225232 85684835928344806933 OD-02 
3. 58221 7514683983264751 53151420-02 
3.60588704656712232295098642910-02 
3. 59760062851346340978567545500-02 
3- 56104273245619168958099261460-02 
3.49960379221006015283231269430-02 
3.41643891193175350272774468440-02 
3. 31450569643689019454045651590-02 
3.19653831532927540866569041820-02 
3.06531230468658719484124912940-02 
2. 923153098624495474421 08882480-02 
2.77244033861403325752035526730-02 
2.61 535939874822903169080855780-02 
2.4539511 5965532382553594368420-02 
2.2901 10786389161583875 15422320-02 
2. 12553607073190895213605882290-02 
1.96197575775782409144376403000-02 
1.8007281 7244301 1026085614396 70-02 
1.643 1403 830946 73646 5 31 2597 11 90- 02 
1.49035807688038169147992147980-02 
1. 3433 762 737352365 98202 11 1326 10-02 
1.20304096501063649487343700190-02 
1. 07005173007415836473956290200-02 
9.44965356067182891410227548200-03 
8.28200462055587478896068258960-03 
7.20043108085676398163356089720-03 
6.20 65 335 163232 368 1 04 1 1 479405 30-03 
5. 30072698211976176990669316660-03 
4.48232379256988813586139105190-03 
3. 74962 3785 1732 40983012 5223 43 30-03 
3. 10001118134606174438582477130-03 
2. 5300570702425291 8197550867010-03 
2.03562647232424776049218663160-03 
1.61198888451271007900590919090-03 
1.25393116834045307589050315580-03 
9.55871616119136634409414321420-04 
7. 1197401 7389998913679991701990-04 
5. 16260540425117749359003586410-04 
3. 627223208858656369171 12772060-04 
2.45426451451307087084098595230-04 
1. 58618551782095418120585579250-04 
9.68195870013249997579195468270-05 
5. 49 16 1 1 0304870 900 19603 15471 8 00-05 
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TABLE m.- Continued 


5? 9.773913170480593^5450292307810-01 2.82429506467925583803161063670-05 

53 9.84877420953869501744323069730-01 1.26575089035486858109269137050-05 

54 9.90893844070220765232994487800-01 4. 60490466968268215759252906680-06 

55 9.95422381760313970000208480680-01 1. 17 33 100 15625722 184585422421 00-06 

56 9.98449732520253187795409873870-01 1.38907604030270876502485429920-07 


N = 64 


M X( M) 

1 2. P42 734756368509 951 43 17 99 17929-04 

2 1.62875564627137524880239805910-03 

3 4. 11023335424740037704506689190-03 

4 7.7 33 2 7492 78 7 0981 796 6903 2 4 66 710-0 3 

5 1 .24842564284341 23755081 5023793-02 

6 1. 836 36 74565 75019 4842 733 4992420-02 

7 2.58667440461597347375406007483-02 

8 3.34486823953361089424180004540-02 

9 4.26218765572054122512535324983-02 

10 5.2 85 599 3 83918763 516743 3 44972 43-02 

11 6.4123081237164263739971 1736240-02 

12 7.641264765893224949 11535074033-02 

13 8.96817456 05 5 5 773 582 99 51 4 862 070-02 

14 I .03905052066138250230559137630-01 

15 1. 19049950963193038280375379640-01 

16 I .3 533 16 18 20 66 36A 980 776 1 36 883 70-01 

17 1 .519631099618593 4056436 9194120-01 

18 1.6965545446091 9455958731009010-01 

19 1.83117 74750895926410888 199041 )-01 

20 2. 0730725172026826457811 0295600-01 

21 2.271 79499436933618168807200900-01 

22 2 .476883992247554083033201 84480-01 

2 3 2.6 87 8634*^ 79813367034631 1510090-01 

24 2 .90424333192796462817546367350-01 

25 3. 12552070979817907005564561680-01 

26 3.35118103294198705338060873630-01 

27 3 .58069930434592060739168971410-01 

2 8 3.813 54132 77 7 55 3103 03 9905 6 50080-01 

2^ 4.049 16496738436533392254896900-01 

30 4. 28702142501904863688808629690-01 

31 4.52655653237412791762170355490-01 

32 4.76721205508834938067500913010-01 

33 5.00842700582442288586044275660-01 

34 5.2496389633359 1981243848494560-01 

35 5.49028539449689722400694128670-01 

36 5.7298049762514914971851 7795070-01 

37 5 .96763891443157682955107082020-01 

30 6. 20323225639024131110360210110-01 

39 6.43603519440697264811452018100-01 

40 6.66550435683679604289153827420-01 

41 6.89110408399993260022540004770-01 

42 7.11230768584063524221191473320-01 

43 7. 32859867842352426060867653580-01 

44 7.53947199638275497046479741450-01 

45 7.74443517849364572376900137690-01 

46 7.94300952359759199857376253920-01 

47 8. 13473121417927811890907316330-01 

48 8. 319 15 24049 70013 330 02 96^4 01050-01 

49 8.495842274032457375 19401678140-01 

50 8.66438803387006843908118631310-01 

51 8. 8? 3 95900 19 895 10778 194 7638 1 60-0 1 

52 8.97549201612197228559911554930-01 

53 9. 11732332955391692676318464110-01 

54 9.24955842186686826770258631540-01 

55 9.37188828586012868265585906410-01 


W(M) 

6.368026980122339277523403802 10-03 
1.22645019723311489026383415100-02 
1.6729^518087043355905445725860-02 
2.02863927661332712943210565230-02 
2. 31661137595694234443737143920-02 
2.5499 5 86 037494242 22 15192 58 8460-02 
2. 73726729956179404216862952310-02 
2.88470969076353420322571518060-02 
2.99703195497889752824717683450-02 
3. 07807677398731637964269838030-02 
3. 13108349946321272172266756230-02 
3. 15 88 7108 7 730284086875645644 80-02 
3. 16395449499336248070299890390-02 
3. 14862107185876677406895400990-02 
3. 11498178928057178165291000190-02 
3. 06500603908746282906371 310830-02 
3. 00054539762574484415684978710-02 
2.92334980206947549678059936960-02 
2. 83507842158630483780882894390-02 
2. 73730677754200720159246552870-02 
2.63 153119864947370866693106820-02 
2. 5191 7138731594526818552559220-02 
2.40157166347653136791420658160-02 
2.28000130648664546408680870920-02 
2. 1556543 122855 7770776500 509430-02 
2.02964880817612943368468562620-02 
1.90302631222219956098970780690-02 
1. 77675098253027863737675236150-02 
1. 65170896957933351703542 7434 10-02 
1.52870795958734977389764352080-02 
1.408476976778264090309763971 50-02 
1.29166649602832507204409916260-02 
1. 1 7 88 4890379749724 1 34460 82 62 60- 02 
1.07051933382587144057193537800-02 
9. 6 7 09 689 431 92 68 76 008930961 98 70-03 
8.68926294915834710973386464500-03 
7. 762798743638852516 71115126770-03 
6.89360023365815431786362022780-03 
6. 08301991317956920806839849400-03 
5.33177060602252863204134786960-03 
4.63996067589223966819026182190-03 
4. 00713245538966920352762 728750-03 
3.43230361098042949054831182030-03 
2.91401 113055570129107292455190-03 
2. 45035759419534124924603248740-03 
2.0390593668742 1758243426689420-03 
1. 6774 963 340048483 4 8205 663 89660-03 
1. 36276278677298663542597866440-03 
1.09171905410178308582268597980-03 
8.61043471684468173120426092610-04 
6.67284275765385991523889695520-04 
5.0691 1010130937422282772 758150-04 
3. 76365038994560573338505211030-04 
2.72100766012640803925741426050-04 
1.90673170429182532403086423050-04 
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TABLE in, - Continued 


56 9.^8^02705132273150632594618970-01 

57 9.58571265668052489631096141920-01 

58 9.67670746559112230147517900650-01 

59 9.75679882812913767643892966720-01 

60 9.82579958827055241937253557810-01 

61 9. 8835485462 1831220627568329040-01 

62 9.92991091279751232305627621980-01 

63 9.964778959333901 55286397207660-01 

64 9.98807475070240595140865748180-01 


1.28703285181654322752992051820-04 
8.30012585609039096913176826040-05 
5.05666605432121617432357641500-05 
2.86337170277638388481023376270-05 
1.47051797303606779765185540620-05 
6.58256620421941589465720261530-06 
2. 3925331 7169766559592308428910-06 
6. 09176852250681144009880866140-07 
7.20862411021308577614486841420-08 


N = 72 


M X( 

1 2.2655850171B84D8 14961236602130-04 

2 1.295945836422682 18060657798560-03 

3 3.26332184588186790921547803370-03 

4 6. 1456 2? 88 274 173 11 5444 7524 58430-03 

5 9.92535940422962514336885081070-03 

6 1.46024350090265862439289550470-02 

7 2.01696095223247435889461686630-02 

8 2.66177018692982062092119672500-02 

9 3. 3935703931 3054750775997159150-02 

10 ^.21108567508967829283140473760-02 

11 5. 1128709056365432681 8654065580-02 

12 6. 09731671787292405863614915180-02 

1 3 7. 16265409400830932244548118780-02 

14 8.30695880037562807393392289370-02 

1 5 9.52815 581 1527223 3 53 840842 88780-02 

16 1.08240238084269502590099350530-01 

17 1.219219980231160 10 14892 5074570-0L 

18 I •36301839155050713522945290660-01 

19 1.51353443376233212320981498690-01 

20 1. 670492246 7 16082 53 59P 85 0435 5 70-01 

21 1. 833 603 8242698 40 72873 645436320-01 

22 2.00256956639300345276445971110-01 

23 2. 177 0788499 77879 5144914 8243 010-01 

2 4 2.356 8 1061786926 2 I IS I 32 75 5 90240-0 1 

25 2.54143398545451273846026866010-01 

26 2. 73 06 08864 026457 461 6 4015 736930-01 

2 7 2.92-3 9 866000 27 383 5342159 97 86 000-01 

28 3.12121062919783202689799591450-01 

29 3.32191714458403725621303945230-01 

30 3. 52573577729950314022619725220-01 

31 3.7322 9028 88 8702 3 9509923 12 38 740-01 

32 3 .941 19927408579980784772198120-01 

33 4. 152 0 768727 7296 8 93 731 78 7931620-01 

34 4.36453348981900292525269934540-01 

35 4.5 7817652 15 7 1 339 67D I 543 8564150-01 

36 4.79261108765987707166130085380-01 

37 5. 007440766801152 79575 701368910-01 

38 5.22226833526491229390768691170-01 

39 5.43669650665710241530791595070-01 

40 5.65032867166695166769852086920-01 

41 5.86276963642257477234596897370-01 

42 S. 07 36 2635 80 993 5 7909 1 862 7 8837 30-0 I 

43 6.28250867642814941355660808300-01 

44 5. 8903 003 97559 19 192 4145 02426 20-01 

45 6.69280822431998853479379086200-01 

46 6. 89346604540810639702181741210-01 

47 7.09063205909049891269663553390-01 

48 7. 28394125322648863513400945370-01 

49 7.47303572646731311 140467064390-01 

50 7.65756535399831375826892676040-01 

51 7 . 837 1 8843 878 765 1 01 33273 7466 1 00-0 1 

52 B.Oll 5723471 3517634029745524360-01 


WIH) 

5. 2102 522 119183 55775351 8556876D-03 
1.00988872373720195635071 399380-02 
1. 38567033625125639387410013300-02 
1.69000329163533301385554513520-02 
1.94119346522640019293614073290-02 
2. 1496 1205193572764210000261050-02 
2.32200736369783365324794464850-02 
2. 46318352109771274514059672790-02 
2.57678986707064465400783080840-02 
2.66574008516322690680086068430-02 
2. 73245429551757731692564772260-02 
2. 77900775599680088893504526380-02 
2. 80722649337077032133382 955170-02 
2. 81875097946625456297244978870-02 
2. 8150796522406917976224401 1010-02 
2.79759923472178831877038331120-02 
2. 76 76061 3334573389676692599160-02 
?. 7263216527127410665250251 5070-02 
2.674902833033766515785745643 lD-02 
2.614^5014314372693190182337360-02 
2. 54601286694612900354335033870-02 
2.4705 928 197778115 90854918 I 7730-02 
2.38914680907009480517975854420-02 
2. 30258818312600638473312496780-02 
2.211 78770061865853237220547390-02 
2. 11757396541996820276290419170-02 
2. 02073340401832289135666765310-02 
1.92201017911331388648256216610-02 
1. 82210584904624027452744146710-02 
1.72167901235100252380006780410-02 
1.621344930377813578588 12042280-02 
1. 5216 752 37587314322000875738 70-02 
1.42319760797721832351632039300-02 
1. 32639568664186487952322326430-02 
1.23170902928029555066551096310-02 
1. 13953323526967255953024561310-02 
1.05022023149288514305075369990-02 
9.64078724291904771799579707780-03 
8. 81374825591961106418261784480-03 
8.02332856315722511699816963930-03 
7.27136327613258551842127671160-03 
6. 55929098121036991429536822420-03 
5. 88816703392393266276282953500-03 
5. 2586 785 1783030067434386436050-03 
4. 67116079405499692844626708360-03 
4. 12 56 1555269043266487628090130-03 
3. 62 17302638245671 09556490003 10-03 
3. 15889891407777651986072739290-03 
2.73624390403826058529177481990-03 
2.35263897285007137707935907060-03 
2.00673300837591306031453131870-03 
1.69697459479122492030331862520-03 
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TABLE in, - Continued 


53 8. l8D394127343553387254^t6779l5O-0l 

54 8.34334111036532053005185140543-01 

55 8. 500 1 1 149 13 121838 532 313 461 8500-01 

56 8.65041489074891711274692528110-01 

57 8.79397289473172661436264844860-01 

58 8.93051957259107471979258535410-01 

59 9.05980197150146930403579955630-01 

60 9. 181 58058692611989133913292900-01 

61 9. 2956 2980807 3S04872 59642988523-01 

6 2 9,401 7 3 833 7S 5 7714 249 92 695 55 006 0-01 

63 9.499 70958 45192 53 33 80066 6454973-01 

64 9. 58936203054805800619407789070-01 

65 9.67052956783109471320429835360-01 

66 9. 743 06 1 809592 89 5 703 5 L 40 5 95752 0-0 1 

67 9.805324372064317 18315052979350-01 

68 9.86169913160150991351358214170-01 

69 9. 90 7 5 8446 2441 739 26 89325 80 18800-0 I 

70 9. 9443 9 54860 7452 4 7407965 8 047340-01 

71 9.9720644935270334037975891 7070-01 

72 9.990 54303 7749099 42 72 965 7010710-01 


1.42163714413491195676022805160-03 
1. 17884445421753928728978311750-03 
9.66596532344110827064620274960-04 
7.82795522520689057577613677790-04 
6.25271573156678723277208375700-04 
4.91808482720940341237697668940-04 
3.80168962331445711777036814030-04 
2. 88119356823933039044610433420-04 
2. 13453669421153828438144593910-04 
1.5401 6739673 94791 14902 73006750-04 
1. 077264298284182722218304802 10-04 
7.2594681 1464571107408430587120-05 
4.67473089230037426933609926840-05 
2.8442412952591 1194038140740090-05 
1.60872889466509307234993183190-05 
8. 2537358353542231 1547211643680-06 
3. 691663944444199272405643561 80-06 
1. 3409196580877728429915421 7790-06 
3. 41253309169770412641426506810-07 
4. 03687308200687593243976524850-08 


N = 80 


M X('^) 

1 I .34988098941977736467724275370-04 

2 1 .055796217729079 31 2476793 642 80-03 

3 2.661 1 9381081 0903 78 I 3 8 1 63 76 9?60-0 3 

4 5 .00299915496987063564552676603-03 

5 8. 0 79 85 63 440 54 2 84 032 349 80 81 246r>-03 

6 I . 183 858 I 399 629 38 491 3 698 4906 540 -02 

7 1. 642451922803655700131161 84890-02 

8 2 . 168 16970 8 3040 76 75 1 42 72 6 06 1 5 70-02 

9 7652908?38l 988870 16 3139324570-02 

10 3.4329765801 1 0333 0 754825 0 5 54890-02 

1 I 4. 170^742161 58419 6864 6687882910-02 

12 4. 976 1201 1090417544024206434060-02 

13 5.849 34261 74637 77 483 7 1730 72 7540-02 

14 6. 7 886647140 1 8 362 88484 9Q? 003 070-02 

15 7.792 706583 70 7372 815 7 3274 127 330-02 

16 8. 869988 IQ 13 33475 34331 157399660-02 

17 9. 98893189904369480533481271 150-02 

18 I .11778651476604537340633447900-01 

19 I .24250232205850937492050497110-01 

20 1 .3723552100831 3360973668300260-01 

21 1.508 651 142745 813 3 700267 2 45 3 520-01 

22 I .64963775546758672801884402580-01 

2 3 1.795754671479773 34736595618040-01 

24 1.945 6 3 382 846 1226 3315 82 2 858 0220-01 

2 5 2. 102 0 9981 536-0808 0 33 8148 2 0906 00-01 

26 2.2619202191 56675001 17676018880-01 

27 2.42585598123358912585737198130-01 

2 8 2.59 3 66 176 30 0275 3 73008 773 518910-01 

29 2 .7650363^0342481 74792135086000-01 

30 2. 93 9 8 72 836 3 99 75 S 42 2 45 8 7 396 7690-0 I 

31 3. 11 7759 56 22 60648 4414345 51 PI 760-01 

32 3.29847971497031902654384005260-01 

3? 3 .431 76238236003075928883218000-01 

34 3. 6h7 3 3 26 84 118520 3666474 6655 930-01 

35 3. 35491223852775152889409582550-01 

36 4.0442 1958426803049989200676470-01 

3 7 4. 2 34 9 706065 9494 5 9 1C9 29 73 3 74460-0 1 

38 4.425878967898118 22964340813770-01 

39 4.61965654012514652078088606420-01 

40 4. 81301384153279089342938416580-01 

41 5.00666047432382723749551218670-01 


W(M) 

4. 3495652060892343014475096 7820-03 
8. 47578195241485945338194644750-03 
1. 1686 107 756646 4 2 74 48 7356 5 36700-02 
1. 43 19415 8109843 2315 5 559 5 5546 60-02 
1.65 24 69 5 163452 38619 5 8904 799410-02 

1. 83 86 174 933222 4 755 8167 52 4003 80-02 
1.99584239200327903973049942100-02 

2. 12 80 05 9409091643 7595813 0013 90-02 
2.2380142 8450629034998596577540-02 
2. 32816601926778140790640447310-02 
2. 4003485 102954722 87 1 08 9757 82 50-02 
2.45616107663741673430094727070-02 
2.49699438666483405882481826230-02 
2. 52408376398618408996199009180-02 
2. 53 8546 OL 76898 72 76 83 55 36 38 9940- 02 
2. 5414054616 388 6402647522 75 496 0-02 
2.53 36 126 102 843 414876859166 76 30-02 
2. 51605777923909806629886178900-02 
2.48958106046957076925960370230-02 
2.45497966905653661668731340710-02 
2. 41 3013354524399592059529631 70-02 
2.36440836911752615862751473550-02 
2. 30986034981585845074822315280-02 
2.25003637728439300487551415790-02 
2.18557640908462901554549203360-02 
2. 11 709423730110247286086082210-02 
2.04517808639311454534359817770-02 
1.97039094169597357456128321620-02 
1.8932 706 7995 83 452 966907993 68 10-02 
1.814330058826464547955422375 10-02 
1. 734056611030631036343267611 10-02 
1.65291248031985488203453513780-02 
1. 57133422930249256160185828160-02 
1.4897326438335337980509911 7440-02 
1.40849255411399426985286762740-02 
1.32797268899230414543208028920-02 
1.24850557690149200902652708140-02 
1.17039750429430995 145963426750-02 
1.09392854024927043814642 797690-02 
1.01935263403637921607336683520-02 
9.46897790792912403904446014350-03 
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TABLE in. - Continued 


4? 5.200 3 0556 3 36 I 576 8 5P3 7582 725 1 60-01 

43 5.39365819643078386864798098180-01 

44 5.58642786373098627256138095720-01 

45 5. 77832489699726051687168939570-01 

46 5.969060907243128 38 23 38 79 367 770-01 

47 6.1 58 349 22 052 146 9 5b0 8 206 04 745 80-01 

48 6.34590531 105145653615305058600-01 

49 6.53144723106279696397544893060-01 

50 6. 71469603671291088545781 557450-01 

51 6.89537620943806541567302131330-01 

5? 7.0732 1607210591842329800706520-01 

53 7.24794819934436594310017737680-01 

54 7.41930982143002673411263691280-01 

55 7.58 7 0 4322 11 291 12 2 50 3325 32 88490-0 I 

56 7 .750896122893802 39389488976130-01 

57 7.91062207382855724519959920380-01 

58 8.06598081585302330408246350260-01 

59 8.21673864850132055560956517070-01 

60 8. 36266878 18 11 505 758 87550 1 06640-01 

61 8.503551678773881 55376952315800-01 

62 8.6391 753863225479502211 7125560-01 

63 8.7693358578711674379720^496170-01 

64 8.893 83 72564 1134578^82612 54863^-01 

65 9.0124922552103512795471 7762790-01 

66 9 . 125 122 318164187 333 38 26 5 647650-01 

6 7 9.2 31 5 5 796 98 82471 769 3 09 5 3 08 3 790-0 I 

68 9. 331639 05 16 03 52 2 6469 007 00956 70 -0 1 

69 9 . 4252 I 49630622 5 T 0735959 2 5 88000-01 

70 9. 5 12 I 44889Q6 1578 396 1540 0463 310-01 

71 9. 5922980167334199648225 1366430-01 

72 9.665 5 5 37243 26828 819 29 76 1263 310-01 

73 9 . 731801 77234825759221 353553300-01 

74 9. 79094246905973605130670543830-01 

75 9.34288681 91 0 1 8 56 5 846 2 1 5 39 1 2260-0 1 

76 9. 8375566687171003170821 8226080-01 

77 9.9248848353464427253516 1932260-01 

78 9.95481525936428309057382422970-01 

79 9.977303303094831 57292129891800-01 

80 9 .99231 741 55 5 5 80 3 0 574063 6 3 34 1 00-0 L 


8.76766329021381732915419525820-03 
8.09135222348729698984011612910-03 
7.44156526850583272030272982640-03 
6.8195 789422592 72536050529925 00-0 3 
6. 22643170189369069071330864400-03 
5.66293076617417343517951 796700-03 
5. 129659752318167850531 60419050-03 
4.62698709919223663084586673220-03 
4. 15507524166065959085427195080-03 
3. 71 389049525709037214282525340-03 
3.30321360528109238817457542170-03 
2.92265090987122223397065464840-03 
2.57164606255168465229899898290-03 
2.24949225617193025798686447490-03 
1.95534488704288105453077960020-03 
1.68823459540737513184857457120-03 
1. 4470806 1615543045608121 792930-03 
1.23070437189781947622930991800-03 
1.03784323913577895672867711190-03 
8.6716441 7302386882028608532490-04 
7.17278829894207742627634624350-04 
5. 86754986752020082087432003420-04 
4. 74 13 2 73 6778163 459 796548 8554 50-04 
3. 77936840936066299789533397050-04 
2. 96690296754927234461100658160-04 
2.28 92 734 549860485106464 7 720840-04 
1. 73206097640938731938015550650-04 
1.2812071 0807468536470570231440-04 
9. 23130594823897408052068092210-05 
6.44838370251631022096360566650-05 
4. 34030338781 156887277212731000-05 
2. 79197390083197977796094495930-05 
1.69712150967909825511176581890-05 
9. 59120169751647037720725821080-06 
4.9174045 191215538 1367442297150-06 
2. 1981 331 13579720262080222635 10-06 
7.980535697742221504 19927165150-07 
2.03027290066793893914701039510-07 
2. 40 116243686377975818236737870-08 


N = 83 


M X( M) 

1 1. 53869471730449127787526061 090-04 

2 8.768 16829827736341633302732980-04 

3 2.20897 156136370891929686091620-03 

4 4. 152 02 ^ 772 5 78924 573 3464 5 748940-03 

5 6. 705? 5 1186507046 14441225542490-03 

6 9. 866 59911 13 316941440436 3 4655 90-03 

7 I .363 2 9 632 17942200439 1 352 5 76770-02 
I .8000 3 0949343207 951 28996 224390-02 

9 2. 29637383329471672481 719637300-02 

10 2.85175245776011167132773130020-02 

11 3. 46551458112921470183431856500-02 

12 4.13693045025970786363872289150-02 

13 4. 86519468 15 788 03 44526 74 2 5482 20-02 

14 5.64942796359466958213486498860-02 

15 5.4836786705399197688531 7271670-02 

16 7.38192444169743612958833281520-02 

17 8.32807 17607368420943971 8718030-02 

13 9. 32 59675572 12162 20585702 26967D-02 

19 1.037438084473059 19103124390330-01 

20 1 . 14720244053354216403624958610-01 

21 1. 261 754652631 76960440139519650-01 

22 1.38095347933800047391453402070-01 


W(M) 

3.69 114735498807923620372742790-03 
7.22563894337416823454146543860-03 
1.0003 151 0794536814061781041530-02 
1.230487891788869587992 19753160-02 
1.42544844072871169083622477000-02 
1.59217925543483316019796944030-02 
1. 73519216240375772605342380620-02 
1.85766479510714579903301389650-02 
1.96197653235890735123358553810-02 
2.0499941 7148971596285303492860-02 
2.12323785800910139383286365850-02 
2. 18298377835222770342850741070-02 
2.23033037760054282483874115710-02 
2. 26624588378135610442970589570-02 
2.29159452268487300341960983490-02 
2. 30716432317583599671658248930-02 
2. 31367934466881118491823437900-02 
2.31181557229325335031948109200-02 
2. 302204375048862599310 13650320-02 
2.28 54442 7915 970676 94170288 3860-02 
2.26210353613729205047014845320-02 
2. 23272445825442512084832937770-02 
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TABLE m. - Continued 


23 1 . 50465179424444811018103619960-01 

24 1.63269677712244121537249388850-01 

25 I .76493011108019308154191034020-01 

26 I .9011 8818567924312745540^83550-01 

27 2.04130230591632933322473757640-01 

28 2. 18509890693709178074727121160-01 

29 2. 33239977431535721281866356600-01 

30 2.48302226970543514655488423870-01 

31 2.63677956165268732647083772900-01 

32 2. 79348086 132922331325235208360-01 

33 2.952931662945791 046601 65289030-01 

34 3. 11493398857749438929883702300-01 

3 5 3.279 2 8663 7129083 780 5033 78 52 960-01 

36 3. 44578543715555824477073919330-01 

37 3.614223503244701 83648958560600-01 

38 3. 7B439 14956 604089764256 2 535380-01 

39 3.95607788293873460104421663370-01 

40 4.129069207122573161 87774867660-01 

41 4.3031 5035 131560372222400073960-01 

42 4.47810480923158873678851914370-01 

43 4.653 71495641122091013 504 881860-0 1 

44 4.829 7 62 322 7 7 5441 938942 I 8080100-01 

45 5. 00602786618147474470302325800-01 

46 5. 18229224664569165009239476610-01 

47 5. 353 33610089 5863 6647859 71 47930-01 

48 5.533940316914283 1338081 3 56496D-01 

49 5. 7088863081 3270883253239893470-01 

50 5.8829562869400364111411 1345310-01 

51 6. 05593353716403760086625659000-01 

52 6.22760268518942807499824047330-01 

53 6.39774996937590880404154651 020-01 

54 6. 56616350744166998874271950540-01 

55 6 .73263356 148014312269987900170-01 

56 6.89695280028052780974279539430-01 

57 7.05891655862579845000465020900-01 

58 7,21832309324550462011336699750-01 

59 7. 374 9738 3510 47 10 6964478 85274 1 0-01 

60 7. 5286 7363 77 14867 7174854 1921590-01 

61 7.6792 3 10211 57 266 2 5236532 261820-01 

62 7. 82645841 151497561630323419870-01 

63 7.970 1723754148245216 1789218850-01 

6^ a. 11019384930701349358425429110-01 

65 3 . 24634836375636247998356253460-01 

66 8 .37846626078593064254873323560-01 

67 8. 50638290680129248713656167820-01 

68 8 . 62 993 88980 298 59 5026 71 48370030-01 

69 8 . 74398025990029557700194876790-01 

70 8. 863358639551 7874872921 3265340-01 

71 8.97293149 131430502369304878580-01 

72 9. 07756 225 49 28536 741 22 174 122 960-01 

73 9. 1 77 I 205262 83956 2531161 2013 530-01 

74 9.27148222046328992142132627690-01 

75 9.36052972689217583880964465140-01 

76 9.44415205640429309424232407570-01 

7 7 9 .522 2 449 8004 5323 5460 3654 8933 50-01 

78 9.59471115945514794536913116500-01 

79 9.66146026868960918216036081130-01 

30 9.722409107377^2600409556164810-01 

81 9. 77 74817051 76633 5U 74 79 7005540-01 

8? 9.825 60 94176 2238 1 1542 471 64 79600-01 

8 3 9. 86 97 3 101384314822741256483340-01 

84 9.906 792 75 762 9882 5935 891 8 420 360-01 

8 5 9 .937 74 848695501760080236 172070-01 

86 9.96255971226792933206866259110-01 

87 9.931 1958 42 5077 73 8 80658? 3564450-01 

88 9.99363554641801509166582086000-01 


2. 19782632140424635420960764430-02 
2. 15 790 75 5252 30 770 7 329518350640-02 
2. 1 1 34 473 64041 995 807 8920 723 55 7D-02 
2. 06490695868463852838046271470-02 
2. 01273039947629874600639988960-02 
1.95734521883905714691899839790-02 
1.89916282497145753049621045630-02 
1. 83857875181945516432640630730-02 
1. 77597278985375354033002890390-02 
1. 71 1709027805550522973631501 80-02 
1, 64613582998089250727961 047960-02 
1. 57953576938748927994147938080-02 
1. 51237553339359291948844061370-02 
1.44480581579288117677832356280-02 
1, 37716120682069409358750065260-02 
1.30971009074155398883643730740-02 
1.24270455901934542629730447790-02 
1. 176380345723958739408899963 10-02 
1. 11 09567906709727483200251 7670-02 
1.04663683479491269792925087930-02 
9.83607051391317609194761589410-03 
9.22037716104487120218594152380-03 
8.62082917867665502964500600670-03 
8. 03880712405830798081268748350-03 
7.47553319376539967498288606550-03 
6,93207363742195314892920082790-03 
6.40934161530294952902935882750-03 
5. 90810049740856753841448854650-03 
5. 42896759797590555286892048440-03 
4.97241833607807350202628257790-03 
4. 53879080992515606472886288390-03 
4. 12829076970279823746883317970-03 
3. 74099697124297620257761528950-03 
3.37686689050261788958447491750-03 
3.03574277671688115237859541500-03 
2. 71 735802018523446582990359120-03 
2. 42 13 43 8089320 582 83716 72 7492 70- 03 
2. 14723604695281 589447947349670-03 
1. 89448250540659305535456074150-03 
1.6624501 7694150295042990070580-03 
1. 45043280233726201434848848450-03 
1.25765853781577888698191497900-03 
1. 08329773070277643154867949060-03 
9.2647077061 8368990639733102270-04 
7. 86255933029293506855983045000-04 
6. 6169753 I 8072 7 36647 27426 U 21 20-04 
5.51 81 3297403820751226559944100-04 
4. 5560 269 7368 53646 00326 58 12 35 70- 04 
3. 72054416202381517548695214720-04 
3.00154011945880701546481037820-04 
2.38 89 136 745106 3 0129 563912715 50-04 
1. 8726 79 5497102 024313 5668764520-04 
1.4^303883517948989116459450100-04 
1.09044699384757816987981652860-04 
8. 05679112918016525920304031260-05 
5. 79892128037689125240748887300-05 
4. 04 68 37 5 956605 8 5443286 72 43 86 70-05 
2. 72147914347288082475234349690-05 
1. 74926321384684394273902141410-05 
1. 06256185742520834287860724780-05 
6.0013661 7874033460350578762110-06 
3. 07529644604226034623260275240-06 
1. 37409555714543299396425677810-06 
4.987061059686907611 33718458380-07 
1.268394732454236635 19623614460-07 
1. 4998454 1449046451303812958200-08 
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TABLE in.- Continued 


N = 96 


M X( M) 

1 1. 30013265739022624601877325810-04 

2 7.39848863076112106089212268540-04 

3 1 .863078977531031 11761955533640-03 

4 3.5012 3323 2 73 7414 71786117411080-03 

5 5.65395492651948724591175064510-03 

6 8. 3198848600564395960471 0345340-03 

7 1.14968903145014230438587509150-02 

8 1.51821612185390723009092186620-02 

9 1.93722^07894267512488742611920-02 

10 2.40631564988311386007736411580-02 

1 1 2.92 50 2413 48908 04 8922 2 65 93 86 860-02 

12 3.492 83499837238802365038037350-02 

13 4. 10917719229286546502468826240-02 

14 4.77342631695736837297560526370-02 

15 5.48490569204787857678965630660-02 

16 5. 2 42 88 7382 3 61 844 772 50 51 06 85250-02 

17 7.04659320640436569941951 685380-02 

18 7. 6951 95746 18 4530 1 75449 05 5 85240-02 

19 3.787819370431841 35136074415810-02 

20 9 .72354127947400881623983120240-02 

21 1 .07013925773626271391515494170-01 

22 1 .17203593750110152231379162750-01 

2 3 1. 2 7793 83926 8262 2 390 3 965^8001 in-01 

24 1. 3877365^0239^0381677376549710-01 

25 1. 501 31630940965553046990453040-01 

26 1.61855936610324700147069363250-01 

27 1 .73934364092193420479257929830-01 

2 8 1. 85 3 5432 6 07 75 7502 5^^3 765 02 in 30-01 

29 1.99102872403801600555393458660-01 

30 2. 121667039164808553777941243 10-01 

31 2.255321 06668532522536123212930-01 

32 2.391 8536644591 7509751042829960-01 

3 3 2.531 I 1983 60 04 115 840 40 69 3 2 64950-01 

34 2.67297488232868630409688821630-01 

35 2 .817270555455989570003745970 10-01 

36 2. 963 3 5 601 62 9 7 786 7 1 9 705 1 6 5 1 6 1 40-0 1 

37 3.11 2 5779939 31 890 5650 5 79 0 755970-01 

38 3. 263 2 80 94 78 103 52 670 3 2 58 4 453 91 D-01 

39 3.415807232181131 90656542555080-01 

40 3. 569997262641628 10256946517340-01 

41 3. 72 56896 846 5B6 53 5 1333 87 41 85500-01 

42 3. 382721 54388602309194859764880-01 

43 4.0409284581079357196972 1 728780-01 

44 4. 20014479063365942686087897790-01 

4 5 3 60 2 03 8? 4966698 129 8 1 554 1 143 30-0 1 

46 4.52093794056958457520080660390-01 

4 7 4.682178789543697 10220626432940-01 

48 4.843 7 5 7474042022 75 5913 72 8002 80-01 

49 5 .005 50472 42 31574 94 0898 42 20 1 950-0 1 

5 0 5. 1 6 72 5107662 1210 8861 786 4 933 700-0 1 

51 5 .32332705256987493301596887810-01 

5? 5. 490063336789811 50082908666990-01 

53 5 .65079095565904092790890402720-01 

54 5 .8108414551 5736735682725183040-01 

55 5.97004707824033614500417386150-01 

56 6. 12824094146637507036062900540-01 

57 6.28525721069 158533824152894800-01 

58 6.440931275650291 78442290204680-01 

59 6.59509992323701350128105621550-01 

60 6.7476015093095741 68653 72 952 080-0 1 

6 1 5 .8932761288331290861078 1105480-01 

62 7.0469657841 8694993968086190760-01 

63 7. 19351455145760894566238350620-01 

64 7.3 3776 8744544258 3765 0 72 70 08 540-01 


W( M) 

3.17546391163804493382470733480-03 
6.24081613331042321461420165710-03 
8.67000343710662532360228251040-03 
1.0700 111837963 2717 3 0385059 52 90-02 
1. 24354490220367850938647527450-02 
1.39348714890301301247900422640-02 
1. 52363431541461114652190767000-02 
1,63665227715754208821856007890-02 
1.73452887665327854942843945330-02 
1. 8188 1537473852219844420986500-02 
1.39076690095791837148917378960-02 
1.951429548152721 554282316379 10-02 
2.00169710249608058193826593390-02 
2. 04234945581274040436581278250-02 
2.07407943612210650378564098590-02 
2.09 75 120 26 592 74 80 7430960643720-02 
2. 113218417481 12524008245281560-02 
2. 121726^5282155111087803014740-02 
2. 12352850167965578004004476710-02 
2. 11908745190530796670964506070-02 
2. 10884131092760956157108572090-02 
2. 0932067572^522767647676119950-02 
2.07258189104126732625230785930-02 
2. 04734836663302101904452808730-02 
2.01737304325111258030671902900-02 
1.9845 09? 5 75 846 38313006 7815 5980-02 
1, 94759779751347591667039180830-02 
1.90746763875268994609437457020-02 
1. 8644 36V9292 1694461 7786641 8050-02 
1. 81881120556504369730892279410-02 
1. 77088803501753316050831288130-02 
1. 72095283710643380562769918330-02 
1.6692811 7607696211377310132270-02 
1. 61613837849438018978329791020-02 
1. 5617 795 4398401064 73 7960044140-02 
1. 50644952434645683814060641980-02 
1. 45038288070075334246240655150-02 
1. 39330382676616546278335438050-02 
1. 336926165120383951 88140702640-02 
1. 2799532222116820462338251 5740-02 
1.22307778701 188726835508852620-02 
1. 1664 320 5 74421 6532 02 86 964 7 54 70-02 
1. 110337598058392264732 9405956D-02 
1. 0548 053 119267201 8286752198780-02 
1.00003542913657921722441066810-02 
9, 4616751 3974928232429394963430-03 
8.93330492411812455921804705230-03 
8.41642701215010333191716609290-03 
7.91 21 1959714174629807210581280-03 
7. 421356649670991406642 15169740-03 
6. 94500910838298183159552413060-03 
6.43 38 463 1279513615203731395690-03 
6, 03 853 76 790026060 21 57945 870820-03 
5. 60965461731813277734354233220-03 
5. 19 76 726 891826 51294522 549652 70-03 
4. 80297399911971522996831602360-03 
4.42 5849 81606568 09 04 21179 1921 40-03 
4. 0665034170786081 1869198127920-03 
3. 72505314520388318818176724270-03 
3.40153567214701201373705247060-03 
3.09590945536845805598064330200-03 
2.8080583 782673 883470608789 2190-03 
2.53779556125711430955931738710-03 
2.28486733075188731046638125900-03 
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TABLE m,- Concluded 


65 7.47957707690393975971^91722530-01 

66 7.61 8790^^2 0 76 72 85 75819332352100-01 

6? 7.75526396365759013133869758080-01 

68 7. 8888533620490148152 1061471200-01 

6 0 8.019413892002674684324763 73580-01 

70 8.14682359662248664903745244820-01 

71 8.2709338301 7598708186818630960-01 

72 3. 391 6 19398728806 89850785608050-01 

73 B. 50875369714514385373158066490-01 

74 3.62221 384231 0383 9894 7728811870-01 

75 8 . 73188080243599653760970831 350-01 

76 8. 83763952231095784180216916310-01 

77 8.93937904436825077071164300550-01 

78 9.036992625439440 58835869226670-01 

79 9. 130 3 77 8490 74955 85 5 5 47 39 9 46710-01 

80 9.21943673331252315625922297550-01 

81 9. 30407583378 125485271747163300-01 

82 9.384206342034241 3 72 6 35 1 9 09784D-0 1 

83 9.45974417900829056513807988290-01 

84 9. 5306 1008 3515954 75 1 7 78 3 1 87491 0-01 

85 9.59672969568273165295618523210-01 

86 9.65303363525245300270654679810-01 

8 7 9.71445 7 57469 8831 897 783 4 8 326140-0 I 

88 9.76594230710724772153450308060-01 

89 9.812433808840371 77231625971860-01 

90 9.853883297119429 16262777261740-01 

91 9. 89024728294633266210906851600-01 

92 9.92148762066895177233592186300-01 

93 9. 947 5 71 55 85 02^24 51 7491 65 869450-01 

94 9.96847180715663512831187493150-01 

95 9.984166718301628 18761264840650-01 

96 9.99464142 121915521487561355240-01 


2. 048957332380055529659471295 lD-03 
1. 82969077411033436802593390600-03 
1. 62 66 38784424324147933 73451770-03 
1.43932237018789928398829583180-03 
1.267219458465338644852054375 80-03 
1. 10976450618183423948942237650-03 
9.66358161271147818974106513290-04 
8.36369458744565163524967680390-04 
7. 19141034983817035429799533820-04 
6. 13993843493168240277797055180-04 
5. 20 2318 553432 5 771949 49092 86090-04 
4. 37146727600282800441352881680-04 
3.64022423157468902500332645450-04 
3. 0313 976 55 700793 1880999502 05 6D-04 
2.44 78 09127 39034090 2547 25.70 7100-04 
1. 9 723 3 73 3 58996 006199043362 0850-04 
1. 5679 6072252520023046220723280-04 
1.22779884657394318044024368800-04 
9.45152329334311844046800196660-05 
7. 13541235032499283847980551870-05 

5. 2 6 74175 339871772 71916 96877660-0 5 
3.78821054575562572395851281540-05 
2.64170193646198372293834003450-05 
1. 7753 494 902188 30642 503 9665 9810-05 
1. 14044486280676190443802272020-05 

6. 92 3 7746 7659692 560 83 832 597 73 30-06 
3.90874683077589680640291643470-06 
2.0021 7537936856465412422231 760-06 
8. 9431 1196291629263145469296280-07 
3.24490369841928467934190394280-07 
8. 25136247374794464169532212170-08 
9.75 5 7 462 101442 911416 56 8490 2 7 20-09 


ERROR ANALYSIS 

This section is concerned with methods for estimating the accuracy of the computed 
quadrature abscissas and weights. The related subject of estimating the accuracy of 
integrals computed by a Gaussian quadrature will not be considered at all because it is 
adequately covered in the literature. (See refs. 7 and 8, for example.) 

Calculation of quadrature abscissas and weights, as in the sample problem, con- 
sists of several computational tasks, any of which can introduce errors into the results. 
For the sample problem these tasks are: 

1. Calculation of classical Gaussian abscissas and weights. 

2. Reduction of the effect of the nonclassical singularity by using Taylor’s theorem. 

3. Calculation of hn and hn by numerical integration using the abscissas and 
weights computed in step 1 above. 

4. Evaluation of bj^^l = hj^^h^ and c^^l = 

5. Calculation of the desired nonclassical Gaussian abscissas and weights from the 
recursion coefficients of step 4. 
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The classical Gaussian calculation of step 1 is very stable. Comparison of abscissas 
and weights computed by subroutine CGAUSS with tabulated values in reference 7 indi- 
cates that the error introduced by CGAUSS is negligible. Subroutine NGAUSS uses 
essentially the same algorithms so it also introduces a negligible amount of error. This 
means that the abscissas and weights computed by NGAUSS will be accurate provided the 
routine is furnished accurate recursion coefficients. Thus, the only source of detectable 
error are steps 2 and 3 above, the singularity reduction and the simpler Gaussian 
quadrature. 

In a computation of this sort it is useful to distinguish between two kinds of error: 
round-off error and truncation error. Round-off error results from performing sequen- 
tial calculations with finite -length computer words, whereas truncation error results 
from terminating a convergent mathematical process after a finite number of steps. For 
example, if an integral is approximated by a finite sum as in equation (1) the truncation 
error is due to using too low a quadrature order, whereas the round-off error increases 
as the order N increases. In a Gaussian quadrature, where all weights are positive, 
the number of decimal digits lost due to round-off is 

I logjo N (59) 

if individual computer words are rounded, and is 


log 10 N 


(60) 


if Individual computer words are truncated. In expressions (59) and (60), N is the num- 
ber of terms in the quadrature sum, which is the product of the quadrature orders if a 
multidimensional quadrature is performed. 

For the sample problem double precision arithmetic is used, thus computer words 
are truncated and expression (60) is appropriate. The maximum N is 2000 so that less 
than four digits are lost due to round-off. Since this loss is negligible it can be inferred 
that essentially all the error in the computed recursion coefficients is due to truncation 
error in the Gaussian quadrature used to compute hn and hjj (eqs. (38)). 

This error can be estimated by using the prc^ram Itself and varying N, For suf- 
ficiently large N, if the leading digits in bn and Cn agree for two different values of 
N these digits can be assumed to be correct, A more reliable way to estimate accuracy 
is to compare with more accurate recursion coefficients generated by a different method. 
The magnitude of the truncation error of a Gaussian quadrature depends upon the order 
and the function being integrated, but is independent of the weight function. This means 
that for the sample problem any value of a can be used to estimate the error. In par- 
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ticular, when a = 1, the weight function is In (1/x) and for this function values of bn 
and cn for n = 1 to 16 accurate to 30 digits are tabulated in reference 7. Compari- 
son with these values shows that for N = 1400 the computed abscissas and weights bn 
and Cn are accurate to about 17 significant digits and the accuracy decreases very 
slowly as the subscript of the recursion coefficient increases. Subprogram COMPARE 
was used for this comparison, 

Langley Research Center 

National Aeronautics and Space Administration 
Hampton, Va, 23665 
July 1, 1975 
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LAGUERRE ITERATION 

Laguerre's iteration formula is used to compute the zeros of orthogonal polynomials 
from their recursion coefficients. The formula is 

x(k+l) = x(k) "P (Al) 

P’ ±\/(n - l)2p*2 _ „(n _ i)pp" 

where p, p’, and p" are evaluated at x = n is the degree of the polynomial 

p(x), x(k) is the kth approximation to the zero being sovight, and the sign on the radical 
is taken to maximize the absolute value of the denominator. Some properties of Laguerre 
iteration are (ref. 9) : 

1. Convergence to a simple zero is cubic. 

2, The iteration is invariant under the bilinear transformation x = (ax + b)/(cx + d), 

3. If all zeros are real and distinct the real axis is divided up into n abutting 
intervals, each containing a zero, such that any starter within the interval will converge 
to the root in that interval. 

4, If all zeros are real, convergence is monotonic. 

Property 1 above means that convergence is rapid. Property 2 makes it possible to com- 
pute the Laguerre iterate of -oo. If the sign on the radical is taken to be the sign of p' 
at the zero, rather than the sign of p* at the current iterate, then the interval of con- 
vergence for the mth zero Xm extends from Xm-l + 0 to Xm+l - 0. 

To describe completely a Laguerre iteration procedure used in a computer program 
it is necessary to state how the iteration is started, what deflation procedure is used if 
any, and how the iteration is terminated. Let x|^^ be the kth approximation to the mth 
zero counting from the left and let x^^ be the starting value. 

Starting Value for xj^ 

The starting value x^^^ is the Laguerre iterate of -oo except for the classical 
Jacobi -Gauss and Laguerre -Gauss cases where it is the Laguerre iterate of -1 and 0, 
respectively. 
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Deflation 

When equation (Al) is used for m > 1, p(x) is replaced by 

q(x) = ^ (A2) 

^ m-1 

n (x - Xi) 
i=l 

and n is replaced by (n + 1 - m). Note that q(x) is not a polynomial but is a rational 
function with a zero/pole pair at each computed zero. Thus deflation has no effect on 
the accuracy of subsequent zeros. This deflation procedure is absolutely reliable as a 
means of preventing convergence to a previously computed zero, provided that no two 
zeros are closer than the resolving ability of the computer floatii^ point number system. 
To apply equation (Al) substituting q for p as given by equation (A2) note that 

m-1 

al = e! y 1 

q P Zy X - x« 

i=l 


and 


£ Ip'Y I 1 

1 ' W/ P IP/^ (X - X1)2 

starting Value for x^^ 

The starting value x^^ for m > 1 is the Newton iterate of q(x)^(x - Xm_l). 

This number can be computed from q'/q and q"/q, both of which are available because 
they were used to compute the converged value of Xni_i. 

Iteration Termination 

The starting value is always left of the zero and convergence is monotonic. Itera- 
tion is terminated when a negative correction is computed because this indicates the iter- 
ation has reached its noise level. 

Other Programing Considerations 

The zeros Xm are computed in DOUBLE PRECISION in the program. To do this 
it suffices to compute only p(x) in DOUBLE PRECISION. The other terms p'(x), p"(x), 

^ (x - xi)"^, and can all be computed in SINGLE PRECISION. 
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SUBPROGRAM DESCRIPTION 

This appendix contains usage descriptions of the four subprograms mentioned in 
the second paragraph of the section entitled "Program Organization." It also contains 
usage descriptions of three other subprograms, JGAUSS, LGAUSS, and DGAMF that 
could be called independently. All other subprograms are described only by their list- 
ings and associated comments in the section entitled "Sample -Problem Program Listing." 

Subroutine CGAUSS 


Language : FORTRAN 

Purpose : To compute the abscissas and weights for a classical Gaussian quadrature. 


Use : CALL CGAUSS(N,X,W, A, B, ALPHA, BETA) 


N 

X 

w 

A,B 

ALPHA, BETA 


Input parameter containing the order of the quadrature 

One -dimensional output array containing the abscissas 

One -dimensional output containing the weights 

Input parameters containing the integration interval delimiters or 
weight -function scale factors 

Input parameters containing weight -function exponents 


Restrictions : X and W must be dimensioned N or larger in the callii^ program. 
In the Jacobi case B must be greater than A. In the Laguerre case B must be 
positive. In the Her mite case A must be positive. All arguments except N must 
be typed DOUBLE PRECISION, 


Discussion : The subroutine is used to compute abscissas and weights for any of the 
three quadrature formulas below: 



- a)^ f(x) dx » ^ 
m=l 


wm f(xm) 


(Jacobi) 
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(x - f(x) dx » Wm f (xm) (Laguerre) 

m=l 

n 

(ax2+hx)f(jjj dx « ^ Wm f(xm) (Hermite) 

m=l 

The Jacobi case is selected by calling CGAUSS with a,^ greater than -1. The 
Laguerre case is selected by setting )3 to a real number less than -1. The Hermite 
case is selected by setting a to a real number less than -1. 

Gaussian quadrature formulas are exact whenever the integrand function f(x) 
is a polynomial of degree (2n - 1) or less. The weights Wm are all positive so that 
numerical stability is assured. The three classical Gaussian quadratures listed 
above have the additional desirable property that the abscissas and weights are partic- 
ularly easy to compute. 

Examples : Suppose the computation of Xm and Wjn for the quadrature formula was 
desired, then 



The call would be 

CALL CGAUSS(8,X,W,0.D,1.D,0.D,.5D) 

The table below lists the last four parameters for some of the better known Gaussian 
quadrature formulas 


Name 

A 

B 

ALPHA 

BETA 

Legendre 

-1. 

1. 

0. 

0. 

Shifted Legendre 

0. 

1. 

0. 

0. 

Chebyshev (1st) 

-1. 

1. 

-.5 

-.5 

Chebyshev (2d) 

-1. 

1. 

.5 

.5 

Jacobi 

-1. 

1. 

a 


Shifted Jacobi 

0. 

1. 

P-q 

q-1. 

Lj^uerre 

0. 

1. 

0. 

-2. 
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Name 

A 

B 

ALPHA 

BETA 

Hermite 

1. 

0. 

-2. 

-2. 

Scaled Hermite 

.5 

0. 

-2. 

-2. 

Gegenbauer 

-1. 

1. 

-.5 

-.5 

Generalized Laguerre 

0. 

1. 

a 

-2. 


Method : In the Jacobi case the abscissas and weights Xm,Wm in 

b ^ 

J (b - x)“(x - a)^ f(x) dx « ^ Win ^(Xm) 

^ m=l 

are computed from the Xm>Wm associated with 

1 ” 
y (1 - x)“(l + x)^ F(x) dx « ^ Wm^Xm) 

m=l 

using the relations 

xm = |(b + a) + l(b - a)xm 


and 


Wm = 


b - a^ 


a!+j3+l_ 

Wm 


The standard-interval abscissas and weights Xm,Wm are computed by a call to sub- 
routine JGAUSS (see JGAUSS write-up for description of method used). 

In the Laguerre case the Xm,Wm associated with 

n 

(x - a)“e-^^^ f(x) dx « ^ Wm f(xm) 

m=l 
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are computed from the Xm,Wm associated with 

n 

x“e"^rO{)dx~ ^ Wm 
m=l 




using the relations 


xm = 


b 


+ a 


and 


Wm = b'“"^e-’^^m 

The standard-interval and scale abscissas and weights Xm>Wm are computed by a 
call to subroutine LGAUSS. 

In the Hermite case the Xm,Wm associated with 

oo 

y f(x) dx ~ Y, Wm f(xm) 

— OO -4 

m=l 

are computed from the Xm,Wm associated with 


f e~x^f(x)cS~ Y wmf(xm) 

m=l 


using the relations 



2a 


and 


Wm = Wm 
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If n is even, let n = 2v and let be associated with 


CO ^ 

y g(|) d| » ^ Wju, g(|jLi) 


M=1 


Then 


"l^+]LL 


= V^ 


Xf+l-fX — \f^ 


Wt»+l-jLl = Wj,+(i = I t^lLi 


(fX = l(l)i^) 


The IjijWjLi computed by a call to LGAUSS with ALPHA = 0.5, 

If n is odd, let n = 2f + 1 and let be associated with 


C ^l/2e-l g(|) d| « ^ o)|ii g(^jLi) 

0 |Ll = l 


Then 


^t'+l+p " \f^ 

— _ — _ 


Xi^+1 = 0 


Wy+l 


yw 

{2u + 1): 


The ^x’e computed by a call to LGAUSS with ALPHA = 0,5. 
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Accuracy : The program computes Xm and Wm to about 92-bit accuracy. 

Subroutines used : JGAUSS, JRECUR, JWTCON, DGAMF, DSQRT, DEXP, DLOG, 
LGAUSS, LRECUR, LWTCON 
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Subroutine PNDER 


Language : FORTRAN 

Purpose : To compute derivatives of orthogonal polynomials defined by three-term recur- 
sion relations. 

Use : CALL PNDER(N,M,X,PD,PS) 

N Degree of the orthogonal polynomial 

M Number of derivatives to be computed starting with the zeroth and going 

through the (M - l)th 

X Argument for which Pn(x) and its derivative are computed 

PD One -dimensional array into which PNDER stores 

lJ--^Pn^^) k = 0,M-l 

K. jjjjK 

PS One -dimensional array into which PNDER stores 

i-^Pn-i(x) for k = 0,M-1 
dx*^ 

In addition to the parameters N, M, and X the calling program must also pass the 
recursion coefficients b^, k = 1,N in labeled COMMON block BOFN and must pass 
ho and c^, k = 2,N in labled COMMON block COFN, 

Restrictions: PD and PS must be dimensioned M or larger in the calling program, 

X, PD, and PS and the contents of BOFN and COFN must be typed DOUBLE 
PRECISION. 

Method : Pj^ and its derivatives are computed from the recursion formulas 


Pq(x) = 0 ■ 

pjj^\x) = 0 


(k = 1,M) 
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Pl(x) = x - bi'] 


p’l(x) = 1 



(k * 2,M) 


Pn(x) = (x - bn) Pn-l(x) - Cn Pn-2(x) 


(n = 2,N) 



(k) 

n 



Cn -L + 1 

y., Pn-2 + (k - 1): 


(n = 2,N; k = 1,M-1) 


Subroutines used: None. 
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DOUBLE PRECISION Function PNFUN 

Language : FORTRAN 

Purpose : To evaluate an orthogonal polynomial from its three-term recursion relation 
Use : F = PNFUN(X,N) 

X Argument from which p^(x) is computed 

N Degree of Pn(x) 

Restrictions : X and PNFUN should be declared DOUBLE PRECISION in the calling 
program. 

Method : The same recursion formula as in PNDER is used. 

Subroutines used: None. 
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Subroutine NGAUSS 


Language: FORTRAN 

Purpose : To compute the abscissas and weights of a Gaussian quadrature formula from 
the recursion coefficients of the associated system of orthogonal pol 5 momials. 

Use : CALL NGAUSS(N,X,W) 

N Number of abscissas and weights to be computed (i.e., the order of the quad- 

rature formula) 

X One-dimensional array into which NGAUSS stores the computed abscissas 

W One -dimensional array into which NGAUSS stores the computed quadrature 

weights 

In addition to the parameter N, the calling program must also pass the recursion 
coefficients bjj, k = 1,N to NGAUSS via a labeled COMMON block called BOFN and 
must pass hp and the recursion coefficients cj^, k = 2,N via a labeled COMMON 
block called COFN. The coefficients bk,ck are used to define the orthogonal poly- 
nomials Pk recursively 

PqW = 1 


pj(x) = x - bj 


Pk(x) = (x - bk) Pk-l(x) - Ck Pk-2(x) (k = 2,N) 


and 



is needed to compute the quadrature weights. 

Restrictions: The arrays X and W must be typed DOUBLE PRECISION and dimen- 
sioned N or larger in the calling program. The contents of labeled COMMON block 
BOFN and COFN must be typed DOUBLE PRECISION. 
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Method: The abscissas are computed as zeros of the polynomial defined by the recur- 
sion formula using Laguerre iteration with partial deflation. The weights are com- 
puted from 


Wk = 


N 

PN-l(Xk) Pk(xk) 


Subroutine used: PNRECUR 
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Subroutine JGAUSS 


Lar^age : FORTRAN 

Purpose: To compute the abscissas and weights for a Jacobi -Gauss quadrature over the 
interval (-1,1). 

Use : CALL JGA*USS(N,X,W, ALPHA, BETA) 

N Input parameter containing the order of the quadrature 

X One -dimensional output array containing the abscissas 

W One -dimensional output array containing the weights 

ALPHA Input parameter containing the exponent of (1 - x) in the weight function 
BETA Input parameter containing the exponent of (1 + x) in the weight function 

Restrictions : X and W must be dimensioned N or larger in the calling program. 
ALPHA and BETA must each be greater than -1. All arguments exact N must be 
typed DOUBLE PRECISION. 

Example : Suppose one wanted to compute Xm,Wm for the quadrature formula 

j 8 

f (1 - x)l/2(i + x)V2 t(x) dx » y Wmt(*m) 

• ,^1 


The call would be 

CALL JGAUSS(8,X,W,.5D,1.5D) 


Method : The abscissas Xm are computed by applying Laguerre iteration to the Jacobi 
polynomial P^“’^\x) to compute its zeros. The function and its first two 

derivatives are computed from equations 22.7.1, 22.8.1, and 22.6.1, respectively, of 
reference (a) of this subroutine. The quadrature weights are computed from 


w„ 
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where 


_ + a)r(n + /3)(2n + a + j3) 

” nir(n + a + 13 - l)(n + Q!)(n + 

Reference : (a) Abramowitz, Milton; and Stegun, Irene A., eds.: Handbook of Mathemati- 
cal Functions With Formulas, Graphs, and Mathematical Tables. NBS Appl. Math. 

Ser. 55, U.S. Dep. Com,, June 1964, 

Subroutines used: JWTCON, JRECUR, DGAMF, DSQRT, SQRT 
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Subroimtine LGAUSS 


Language ; FORTRAN 

Purpose : To compute the abscissas and weights for a Laguerre -Gauss quadrature over 
the interval (0,°°). 

Use : CALL LGAUSS(N,X,W, ALPHA) 

N Input parameter containing the order of the quadrature 

X One -dimensional output array containing the abscissas 

W One -dimensional output array containing the weights 

ALPHA Input parameter containing the exponent of x in the weight function 

Restrictions : X and W must be dimensioned N or larger in the calling program. 
ALPHA must be greater than -1. All parameters except N must be typed DOUBLE 
PRECISION. 


Example: Suppose one wanted to compute Xm,Wni for the quadrature formula 



25 

I 

m=l 


Wm f(xm) 


The call would be 

CALL LGAUSS(12,X,W,.5D) 

Method : The abscissas Xm are computed by applying Laguerre iteration to the 
Laguerre polynomial L^“)(x) to compute its zeros. The function l|,“^(x) and its 
first two derivatives are computed from equations 22.7.12, 22.8.6, and 22.6.15, respec- 
tively, of reference (a) of this subroutine. The quadrature weights are computed from 


Wm = 


xcn 


lS(x) 


n 2 
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where 


_ r(n + a) 


Reference : (a) Abramowitz, Milton; and Stegun, Irene A., eds.: Handbook of Mathemati- 
cal Functions With Formulas, Graphs, and Mathematical Tables. NBS Appl. Math. 

Ser. 55, U.S. Dep. Com,, June 1964, 


Subroutines used: LWTCON, LRECUR, DGAMF, DSQRT, SQRT 
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DOUBLE PRECISION Function DGAMF 


Lai^uage : FORTRAN 

Purpose : To compute a DOUBLE PRECISION gamma function for a real argument. 

Use : GAMX = DGAMF (X) 

X Argument at which F(X) is to be computed 

Restrictions : Both X and DGAMF must be typed DOUBLE PRECISION and X must 
be less than 145.0 and not a negative integer. 

Method : The gamma function satisfies the functional equation 

F(X + 1) = XF(X) 

If X is an integer, F(X) is computed from F(0) = 1 by repeated application of the 
functional equation. If X is a half integer, F(X) is computed from F(l/2) = \/¥. 
Otherwise F(X) is computed from F(Z) where 

(Z-X) is an integer 

144.0 S Z < 145.0 

r(Z) is given by equation 6.1.48 in reference (a) of this subroutine. 

Reference : (a) Abramowitz, Milton; and Stegun, Irene A., eds.: Handbook of Mathemati- 

cal Functions With Formulas, Graphs, and Mathematical Tables. NBS Appl. Math. 

Ser. 55, U.S. Dep. Com., June 1964. 

Subroutines used: DEXP, DLOG 
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